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FOREWORD 


This  report  describes  the  work  performed  by  Systems  Control  Technology, 
Inc.  (SCT)  under  Contract  No.  N00421-81-C-0289  for  the  Strike  Aircraft  Test 
Directorate  of  the  NAVAIRTESTCEN  and  under  Contract  No.  N00014-78-C-0641  for 
the  Office  of  Naval  Research.  The  objective  of  the  initial  work  was  to 
develop  a  system  identification  methodology  suitable  for  extracting  a 
nonlinear  aerodynamic  data  base  from  flight  test  measurements.  The  second 
objective  of  this  work  was  to  apply  this  methodology  to  the  development  of  an 
aerodynamic  data  base  for  an  F-4S  simulator  at  the  NAVAIRTESTCEN. 

Mr.  Roger  Burton  served  as  the  technical  monitor  for  this  work.  The 
program  manager  was  Mr.  James  Vincent.  Mr.  Norm  Franklin  was  principally 
responsible  for  producing  the  system  identification  results.  Report 
preparation  efforts  were  directed  by  Ms.  Clare  Walker,  with  assistance  from 
Mrs.  Toshi  Furukawa. 


"s-  w  .  V.  .  *. 


r 


i 


I  . 


r  «• 


ax»  ay>  az 

bw 

b(-) 

cw 

CO 

CD(.) 

cL 

C£ 

c£(.) 

Cn 

cn(.) 

Cx 

Cx(.) 

cy 

cy(.) 

cz 

CZ(.) 

Fx.  Fy»  F z 

Fax»  Fay»  Faz 
F$x»  Fsy»  Fsz 
fTx»  FTy,  fTz 
f6 
FI 


LIST  OF  SYMBOLS 

Axial,  lateral,  and  vertical  accelerations  ft/sec^ 
Wing  span  (38.67  ft)  ft 

Bias  for  the  (*)  sensor 

Wing  mean  aerodynamic  chord  (16.04  ft)  ft 

Drag  coefficient 

aCo/3(.) 

Lift  coefficient 
Rolling  moment  coefficient 
aCJl/a(  - ) 

Yawing  moment  coefficient 
aCn/a ( - ) 

Axial  force  coefficient 
3Cx/8(.) 

Side  force  coefficient 
3Cy/3(*) 

Normal  force  coefficient 
3Cz/3(-) 


Total  applied  axial,  lateral,  and  vertical  lb 
forces  in  body  axis 

Aerodynamic  forces  in  body  axis  lb 
Aerodynamic  forces  in  stability  axis  lb 
Propulsion  forces  in  body  axis  lb 
Gross  thrust  lb 
Inlet  ram  force  lb 


y 


t 


V 


LIST  OF  SYMBOLS  (Continued) 


£ 


i 

I 

! 


Fn 


9 

h 

hP 

lx.  *y»  ^2 

Ixz 

J 


l  ,1,1 

x’  y*  z 


i 


M 

M 

Mx,  My,  Mz 


mAx»  MAy»MAz 
Msx»  MsytMsz 
mTx»  MjyjMjz 


m 


Net  thrust 

Acceleration  due  to  gravity 
Altitude 

Pressure  altitude 

Moments  of  inertia 

Cross  moment  of  inertia 

Matrix  of  partial  derivatives  of  measure¬ 
ments  with  respect  to  parameters 

Angle  of  attack  numerator  gain 

Noseboom  location  relative  to  c.g. 

Accelerometer  location  relative  to  c.g. 

Net  thrust  moment  arms  relative  to  MRC 
Gross  thrust  moment  arms  relative  to  MRC 
Inlet  ram  force  moment  arms  relative  to  MRC 

Mach  number 

Fisher  information  matrix 
Total  applied  rolling,  pitching,  and  yawing 
Aerodynamic  moments  in  body  axis 
Aerodynamic  moments  in  stability  axis 
Propulsion  moments  in  body  axis 
Aircraft  mass 


ma 

NX 


Engine  air  mass  flow 
Axial  load  factor 


lb 

ft/sec2 

ft 

ft 

slug-ft2 

slug-ft2 


ft 

ft 

ft 

ft 

ft 


lb-ft 

lb-ft 

lb-ft 

lb-ft 

slug 

slub/sec 

g's 


vl 


t 

• 

LIST  OF  SYMBOLS  (Continued) 

Ny 

Lateral  load  factor 

g's 

B 

NZ 

Vertical  load  factor 

g's 

• 

„*• 

p 

Roll  rate 

rad /sec 

Q 

Pitch  rate 

rad /sec 

? 

A 

q 

Dynamic  pressure 

lb/ft2 

?  . , 

R 

Roll  rate 

rad/sec 

Sw 

Wing  surface  (530  ft?) 

ft2 

' 

s 

LaPlace  operator  (o  ±  ju) 

- 

» 

t 

time 

sec 

. 

u 

System  input 

- 

r 

V 

Inertial  velocity  along  y  body  axis 

ft/sec 

9  ..  . 

VT 

True  airspeed 

ft/sec 

Var(  ) 

Covariance  of  (  ) 

- 

i 

W 

Inertial  velocity  along  z  body  axis 

ft/sec 

9,  ... 

w 

Process  noise  vector 

- 

Wf 

Fuel  flow 

lbs/hr 

■ 

X 

System  state 

- 

W  '  '  ' 

t-.-r- 

y 

System  output 

- 

;";-v  5 

GREEK 

a 

Angle  of  attack 

rad 

• 

b 

Angle  of  sideslip 

rad 

Y 

Flight  path  angle 

rad 

— 1 

5 

Control  surface  deflection 

rad 

9—2 

«A 

Effective  aileron  deflection  (spoiler  effec¬ 
tiveness  included) 

rad 

’  - 

r 

L. 

vi  1 

*  ’  *".< 

r 

1  _ 

LIST  OF  SYMBOLS  (Continued) 

«LES 

Leading-edge  slat  deflection 

ra< 

«R 

Rudder  deflection 

rai 

«S 

Stabilator  deflection 

ra< 

«TEF 

Trai ling-edge  flap  deflection 

rai 

Euler  pitch  angle 

rai 

®T 

Thrust  line  inclination  relative  to  x-axis 

rai 

*k 

Measurement  error  of  measurement  K 

- 

p 

Measurement  autocorrelation 

- 

p 

Air  density  lb 

sec2/ft 

°k 

Root  mean  square  error  of  measurement  K 

- 

T 

Time  delay 

msi 

i 

Euler  roll  angle 

ra 

4 > 

Euler  yaw  angle 

ra 

SUBSCRIPTS 

A 

Airspeed  component  in  body  axes 

AA 

Accelerometer  location 

eg 

center  of  gravity  location 

L 

Left 

m 

Measurements 

MRC 

Moment  reference  location 

LIST  OF  SYMBOLS  (Concluded) 


SUBSCRIPTS  (Continued) 

NB  Nose  boom 

R  Right 

REF  Body  reference  location 

x,  y,  z  Axial,  lateral,  and  vertical  body  axes 

W  Wind  component  in  body  axes 


ABBREVIATIONS 

AMB 

Ambient  conditions 

CR 

Cruise 

HSI 

High  stability  index 

FCS 

Flight  control  system 

LSI 

Low  stability  index 

LSI* 

Low  stability  index  (for  comparison  with  HSI) 

PA 

Power  approach 

TO 

Take-off 

TABLE  OF  CONTENTS 


Page 

INTRODUCTION  AND  SUMMARY  .  1 

1.1  Introduction  .  1 

1.1.1  Motivation  for  Nonlinear  System  Identification  .  1 

1.1.2  F-4S  System  Identification  Research  Program  Overview  2 

1.2  Report  Summary  .  5 

SYSTEM  IDENTIFICATION  METHODOLOGY  .  7 

2.1  Overview  of  an  Integrated  System  Identification  Procedure  ..  7 

2.1.1  Flight  Data  Processing  and  Analysis  . 9 

2.1.2  Model  Structure  Determination  .  9 

2.1.3  Parameter  Identification  . 9 

2.1.4  Model  and  Parameter  Verification  .  10 

2.2  Data  Processing  Algorithms  .  10 

2.2.1  Equation  Error  .  10 

2.2.2  Output  Error  . 13 

2.2.3  Combined  State  and  Parameter  Estimation  .  14 

2.3  Some  Practical  Considerations  .  15 

2.3.1  Assumptions  Regarding  Measurement  Noise  Statistics  ..  15 

2.3.2  Number  of  Independent  Parameters  in  the  Model  .  15 

2.3.3  Extrapolation  of  Results  .  16 

2.3.4  Excitation  of  All  System  Modes  .  16 

2.3.5  Effective  Use  of  Sequential  Data  Processing  Schemes  16 

2.3.6  Process  Noise  .  16 

2.3.7  Initialization  .  17 

2.3.8  Numerical  Methods  .  18 

2.4  Considerations  for  Nonlinear  Parameter  Identification  .  21 

2.4.1  Problem  Definition  .  21 

2.4.2  Derivative-Free  Minimization  of  the  Negative  Log 

Likelihood  Function  .  26 

2.4.3  Autocorrelated  Measurement  Errors  .  28 

2.4.4  Spline  Model  Structure  .  31 

2.5  Modeling  Approach  .  35 


xi 


TABLE  OF  CONTENTS  (Continued) 


Page 


III.  FLIGHT  TEST  PROGRAM  DESCRIPTION  .  39 

3.1  Test  Aircraft  Description  .  39 

3.1.1  Aircraft  Configuration  . 41 

3.1.2  Aircraft  Geometry  .  45 

3.2  Test  Instrumentation  .  48 

3.3  Test  Flight  Conditions  .  48 

3.3.1  Phase  I  .  48 

3.3.2  Phase  II  .  52 

3.3.3  Phase  III  .  52 

IV.  FLIGHT  TEST  DATA  ANALYSIS  AND  PREPROCESSING  .  57 

4.1  Flight  Data  Evaluation  and  Error  Correction  .  57 

4.2  Data  Filtering  .  58 

4.3  Data  Reconstruction  .  63 

4.4  Instrument  Calibration  .  64 

V.  SYSTEM  IDENTIFICATION  RESULTS  .  67 

5.1  Longitudinal  Characteristics  and  External  Store  Effects  ....  67 

5.1.1  General  .  67 

5.1.2  Instrumentation  Considerations  .  68 

5.1.3  Results  .  70 

5.2  Longitudinal  Characteristics  at  High  Angle  of  Attack  .  108 

5.3  Lateral-Directional  Characteristics  .  118 

5.3.1  General  . 118 

5.3.2  Instrumentation  Considerations  .  118 

5.3.3  Combined  Initial  and  Final  Estimates  .  119 

5.4  Propulsion  System  Modeling  and  Identification  .  135 

5.5  Model  Validation  . 135 

5.6  Flight  Measurement  Predictions  .  141 

VI.  SUMMARY  AND  CONCLUSIONS  .  147 

APPENDIX  A:  F-4S  MATHEMATICAL  MODEL  .  149 

REFERENCES  .  175 


I.  INTRODUCTION  AND  SUMMARY 


1.1  INTRODUCTION 


1.1.1  Motivation  for  Nonlinear  System  Identification 


'he  Naval  Air  Test  Center  (NATC)  and  Systems  Control  Technology,  Inc. 
(SCT)  have  worked  jointly  to  develop  an  advanced  flight  test  data  processing 
technique  that  supports  an  integrated  flight  testing  procedure  (i.e., 
extraction  of  test  data  for  multiple  test  requirements  from  common  flight 
conditions).  This  data  processing  technique  is  commonly  referred  to  as  system 
(or  parameter)  identification.  iThe  development  of  the  system  identification 
methodology  has  been  pursuedjry  a  number  of  organizations  during  the  last 
decade  [1-7]. _ - 

^Realization  of  this  goal  for  an  integrated  flight  testing  procedure  is 
dependent  on  the  ability  to  identify  nonlinear  aerodynamic  characteristics  and 
propulsion  system  performance  from  flight  test  data.  The  identified  models 
can  be  used  to  define  performance,  stability  and  control,  and  unaugmented 
airframe  dynamic  characteristics  of  the  aircraft  being  evaluated.*  By 
identifying  the  nonlinear  aerodynamic  models  in  a  multivariable,  t®le-look-up 
format,  direct  correlations  can  be  made  with  preflight  aerodynami<f predictions 
(e.g.,  wind  tunnel  data)  and  simulation  models. 
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By  using  a  data  processing  technique  that  can  identify  aerodynamic  and 
installed  propulsion  models  from  many  large-amplitude  dynamic  test  conditions, 
it  is  possible  to  enhance  the  test  productivity  through  a  reduction  in 
required  test  time.  For  the  dynamic  maneuvers,  the  test  time  is  defined  in 
terms  of  seconds,  compared  to  minutes  for  static  tests.  Other  motivating 
factors  that  support  the  development  of  this  technology  include  improvements 
in  safety  of  flight,  and  a  general  expansion  of  requirements  for  higher 
fidelity  aerodynamic  models  of  the  aircraft. 


Safety  of  flight  can  be  enhanced  during  a  flight  envelope  expansion  test 
program  by  using  system  identification  techniques  to  validate  the  aircraft 
simulation  model  for  flight  regimes  already  tested.  The  updated  mathematical 


model  can  then  be  used  to  make  preflight  predictions  for  flight  envelope 
expansion  test  conditions.  In  addition,  when  nonlinear  identification  models 
and  identification  techniques  are  used,  the  pilot's  task  is  greatly  simplified 
since  he  is  not  required  to  maintain  small  perturbation  flight  about  a  trimmed 
operating  point. 


^^The  need  for  improved  modeling  of  aircraft  aerodynamic  characteristic  has 
been,  and  continues  to  be  apparent  in  numerous  areas  of  technical  and 
operational  importance.  Four  such  areas  are:  1)  flying  quality  military 
specification  compliance  testing,  2)  training  simulations,  3)  design  methods 
for  specification  of  aircraft  characteristics,  and  4)  the  development  of 
mission  profiles  that  make  optimum  use  of  the  airplane's  capabilities.  In 
general,  there  is  a  need  for  an  improved  understanding  of  an  airplane's 
aerodynamic  characteristics  to  support  design  improvements  for  increased  cost 
effectiveness,  expanded  mission  flexibility  and  enhanced  operational  safety. 
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1.1.2  F-4S  System  Identification  Research  Program  Overview 


NAVAIRTESTCEN  initiated  its  program  to  develop  advanced  s.  >t 
identification  techniques  in  1971  in  a  joint  program  with  the  Office  of  Naval 
Research,  Naval  Air  Systems  Command  (NAVAIRSYSCOM),  and  Systems  Control,  Inc. 
The  original  purpose  of  this  development  effort  was  not  directed  towards 
enhancing  simulation  fidelity,  but  rather  to  improve  flight  test  stability  and 
control  analysis  techniques.  However,  system  identification  was  quickly 
recognized  as  an  ideal  analysis  tool  for  improving  the  accuracy  of  the 
aerodynamic  data  bases  used  in  aircraft  simulations.  Thus,  during  the  1970's, 
NAVAIRTESTCEN  applications  of  system  identification  were  equally  divided 
between  stability  and  control  analysis,  and  enhancing  the  fidelity  of  aircraft 
simulatons.  These  simulation  improvement  programs  were  motivated  by  either: 

(1)  the  requirement  to  improve  the  aerodynamic  data  bases  of  Navy 
fleet  aircraft  so  that  accurate  simulations  were  available  to 
study  and  solve  operational  problems;  or 
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(2)  the  requirement  to  provide  an  accurate  aerodynamic  data  base  so 
that  Navy  operational  trainer  simulated  flight  characteristics 
would  be  representative  of  actual  fleet  aircraft. 

Development  of  a  nonlinear  system  identification  capability  has  been 
pursued  through  a  multiphase  research  program  spanning  1979  through  1982.  The 
NAVAIRTESTCEN  F-4S  (BUNO  286)  has  been  used  throughout  the  research  program. 
The  test  aircraft  and  its  instrumentation  have  remained  fairly  constant 
throughout  the  research  program.  The  F-4S  parameter  identification  program 
included  three  separate  phases,  which  are  summarized  as  follows. 

Each  phase  had  a  separate  objective  and  its  own  set  of  flight  test  data. 

Phase  I 

Flight  test  data  for  the  Phase  I  research  program  were  generated  during 
the  proof  of  concept  flight  test  program  for  the  F-4S  design.  The  objective 
in  processing  the  Phase  I  data  was  to  assess  the  general  quality  of  the  F-4S 
flight  test  data.  The  purpose  of  this  analysis  was  to  generate  a  basis  for 
formulating  a  flight  test  plan  and  instrumentation  specification  for  the  Phase 
II  program. 

Flight  data  processing  from  this  phase  provided  valuable  insight  for 
subsequent  phases.  Useful  parameter  identification  results  were  not  obtained 
from  the  Phase  I  flight  test  program  because  of  the  developmental  status  of 
the  nonlinear  system  identification  technology. 

Phase  II 

The  purpose  of  the  Phase  II  flight  test  program  was  to  generate  a 
comprehensive  flight  test  data  base  which  could  be  used  to  develop  and 
validate  nonlinear  system  identification  analysis  techniques.  The  Phase  II 
flight  test  plan  was  based  on  experience  gained  from  the  first  program  phase. 

The  flight  test  program  was  conducted  in  a  clean  configuration  with  the 
throttles  fixed  for  each  maneuver  and  with  flight  near  and  beyond  stall.  Test 
conditions  were  generally  initiated  from  wing-level,  constant  altitude  flight 
with  a  ss  10*  and  M  %  .6.  For  some  test  conditions,  the  pilot  applied 
variable  aft  stick  to  maneuver  the  aircraft  through  the  desired 
angle-of-attack  test  range.  For  some  of  the  stall  entry  conditions  lateral 
stick  and/or  pedal  doublets  were  combined  with  the  longitudinal  stick 


command.  Other  flight  conditions  included  single  axis  and  multi-axis 
sequented  doublet  inputs.  These  Inputs  were  made  by  the  pilot  and  were  not 
intended  to  be  repeatable  nor  precise  with  regard  to  their  spectral 
characteristics,  but  were  generally  effective  in  creating  large-amplitude 
motions.  The  overall  test  goal  was  to  force  the  aircraft  through  a  broad 
range  of  test  conditions.  The  level  of  excitation  of  primary  test  variables 
achieved  during  the  F-4S  flight  test  program  is  summarized  as  follows: 

-  Angle  of  attack:  -1*  <  a  <  40* 

-  Sideslip:  |a|  <  18* 

-  Mach  No:  M  <  .6 

-  Rotational  Rates:  |Pj  <  90* /S,  |Q|  <  20* /S, 

| R|  <  25* /S 

-  Full  amplitude  control  inputs 

The  nonlinear  system  identification  data  processing  techniques  were 
validated  through  the  processing  of  the  Phase  II  flight  test  data.  Useful 
aerodynamic  data  were  obtained  from  the  Phase  II  program.  The  lateral 
directional  aerodynamic  characteristics  presented  in  this  report  are  primarily 
from  this  analysis.  The  low-speed,  longitudinal  model  for  the  clean 
configuration  was  also  extracted  from  the  Phase  II  flight  test  data  set. 

Phase  III 

The  objective  of  the  third  program  phase  was  to  generate  an  aerodynamic 
data  base  for  an  F-4S  simulation  at  the  NAVAIRTESTCEN.  This  updated 
simulation  model  has  been  used  to  evaluate  aft  c.g.  flying  qualities  of  the 
F-4S  with  a  high  stability  index  store  configuration.  These  simulation  tests 
were  used  to  evaluate  the  consequences  associated  with  moving  the  aft  c.g. 
limit  further  aft. 

The  Phase  III  F-4S  flight  test  program  included  three  flap  configurations 
(power  approach,  take-off,  and  cruise)  and  two  loading  configuraitons  (clean 
and  a  high  stability  index  store  configuration).  The  flaps-up  flights  were 
flown  for  .7  <  M  <  .9.  All  flight  conditions  were  flown  at  low  to  moderate 


angles  of  attack  (0*  <  a  <  24*).  Sinusoidal  stick  and  pedal  pumping  inputs 
covering  a  range  of  frequencies  were  used  to  excite  the  airplane. 

A  majority  of  the  aerodynamic  models  presented  in  this  report  are  from 
the  analysis  of  the  Phase  III  f^ght  test.  An  important  aspect  of  the  Phase 
III  system  identification  data  processing  task  was  a  demonstration  of  the 
ability  to  produce  flight-validated  aerodynamic  models  in  a  reasonable  period 
of  time. 

1.2  REPORT  SUMMARY 

Section  II  provides  a  technical  background  for  the  presentation  described 
by  this  report.  An  overview  of  an  integrated  system  identification  procedure 
is  presented.  This  is  followed  by  a  discussion  of  several  algorithms  which 
can  be  used  for  system  identification.  Finally,  a  methodology  for 
representing  nonlinear  models  is  presented. 

Section  III  describes  the  test  aircraft,  the  flight  test  instrumentation 
system,  and  the  test  conditions  evaluated  for  the  F-4S  identification  program 
described  by  this  report.  The  descriptions  presented  in  Section  III  are 
provided  as  background  information  for  the  discussions  on  model  structure 
determination  and  parameter  identification. 

Section  IV  describes  the  steps  involved  in  the  analysis  and  preprocessing 
of  flight  test  data.  This  is  a  major  element  of  the  integrated  system 
identification  procedure.  The  overall  objectives  of  this  task  are  to  evaluate 
the  suitability  of  the  flight  test  data  for  identifying  aerodynamic 
characteristics  and  to  develop  a  complete  set  of  kinematically  consistent 
measurements. 

Section  V  presents  the  aerodynamic  and  propulsion  model  data  which  were 
identified  from  both  the  Phase  II  and  Phase  III  flight  test  programs.  This 
data  base,  along  with  the  mathematical  model  presented  in  Appendix  A,  comprise 
a  simulation  of  the  F-4S. 


Conclusions  for  this  study  are  presented  in  Section  VI. 


II.  SYSTEM  IDENTIFICATION  METHODOLOGY 


This  section  provides  a  technical  background  for  the  presentation 
described  by  this  report.  An  overview  of  an  integrated  system  identification 
procedure  is  presented.  This  is  followed  by  a  discussion  of  several 
algorithms  which  can  be  used  for  system  identification.  Finally,  a 
methodology  for  representing  nonlinear  models  is  presented. 

2.1  OVERVIEW  OF  AN  INTEGRATED  SYSTEM  IDENTIFICATION  PROCEDURE 

The  capability  of  processing  large  quantities  of  data  from  flight  tests 
has  led  to  the  parallel  development  of  data  processing  algorithms,  which 
greatly  increase  the  amount  of  useful  knowledge  that  can  be  extracted  from  the 
data.  These  algorithms,  based  on  dynamical  and  statistical  principles,  yield 
very  precise  information  about  the  characteristics  of  the  data  and  the  system 
that  produced  it.  This  methodology  (often  called  system  identification)  has 
been  extensively  documented  in  previous  technical  reports  and  technical 
journals  [1-7].  Key  contributions  in  this  field  include: 

-  The  development  of  an  integrated  system  identification  procedure 
that  includes  a  model  structure  determination  phase. 

-  The  use  of  system  identification  for  supporting  flight  test 
planning:  instrumentation  selection  (type  and  accuracy), 
telemetry  requirements  (sample  rate,  bandwidth,  and  digitization) 
and  input  design. 

-  Algorithm  improvements  for  the  implementation  of  maximum 
likelihood  techniques. 

The  steps  in  the  integrated  system  identification  procedures  are 
illustrated  in  Figure  2.1  with  the  salient  objectives  of  each  step  noted. 

Five  basic  steps  are  shown:  Flight  Data  Processing  and  Analysis,  Model 
Structure  Determination,  Parameter  Identification,  Model  and  Parameter 
Validation,  and  Flight  Test  Planning.  The  following  discussion  reviews  the 
first  three  steps  shown  in  Figure  2.1.  Model  and  parameter  validation 
techniques  are  addressed  later  in  this  report  and  examples  of  the  application 
of  system  identification  to  flight  test  planning  are  presented  in  Refs.  2,  4, 
and  5. 


SYSTEM  MODEL 


Figure  2.1  Integrated  System  Identification  Procedure 


2.1.1  Flight  Data  Processing  and  Analysis 


The  processing  and  analysis  of  flight  test  data  is  a  major  element  of  the 
integrated  system  identification  procedure.  The  overall  objectives  of  this 
task  are  to  review  measurement  excitation,  remove  wildpoints,  reconstruct 
unmeasured  quantities  (i.e.,  acceleration  of  the  aircraft  at  the  center  of 
gravity),  and  develop  a  set  of  kinematically  consistent  measurements.  Where 
measurement  consistency  cannot  be  established,  requirements  for 
instrumentation  system  error  source  modeling  are  defined.  Measurement 
consistency  has  a  significant  impact  on  parameter  identification  accuracy 
since  unaccounted-for  errors  will  bias  parameter  estimates. 

2.1.2  Model  Structure  Determination 


Model  structure  determination  is  the  next  step  of  the  integrated  system 
identification  procedure.  The  purpose  of  this  step  is  to  identify  significant 
terms  of  a  math  model  and  to  generate  an  initial  estimate  of  the  parameters 
for  the  parameter  identification  program.  By  identifying  a  significant  model 
structure  first,  the  potential  for  having  a  divergent  solution  with  the 
parameter  identification  optimization  algorithm  due  to  "overparameterization" 
is  reduced.  Furthermore,  the  cost  of  running  the  parameter  identification 
program  is  reduced  when  good  initial  parameter  estimates  are  used  since  fewer 
program  iterations  are  required  for  convergence. 


2.1.3  Parameter  Identification 

Parameter  identification  is  the  final  data  processing  step  in  the 
integrated  system  identification  procedure.  For  this  step,  refined  parameter 
estimates  are  obtained  by  using  a  maximum  likelihood  algorithm  with  the  model 
structure  identified  from  the  previous  step.  There  are  several  reasons  for 
using  a  maximum  likelihood  algorithm  to  refine  the  parameter  estimates. 
Equation  error  techniques,  which  are  used  to  identify  the  model  structure, 
are  formulated  to  minimize  the  squared  error  difference  between  reconstructed 
and  estimated  equation  variables  (i.e.,  CL  and  Cm,  etc.)  along  the  flight 
trajectory.  If  any  of  the  measurements  are  in  error  (e.g.,  bias  or  scale 
factor),  the  resulting  parameter  values  will  also  be  in  error.  The  maximum 
likelihood  program,  on  the  other  hand,  minimizes  the  squared  error  of  the 
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estimated  measurements  by  estimating  both  aerodynamic  and  instrumentation 
error  parameters. 


2.1.4  Model  and  Parameter  Verification 


The  validity  of  the  parameter  estimates  can  be  established  from  three 
different  considerations. 

(1)  Engineering  judgement:  Are  the  estimates  reasonable  from  the 
point  of  view  of  general  agreement  with  preflight  predictions? 

(2)  Estimation  uncertainty:  What  is  the  magnitude  of  the  +2 a  bands 
about  the  estimate? 

(3)  Prediction  accuracy:  How  well  does  the  identified  model 
predict  flight  test  measurements  for  many  test  maneuvers?  Does 
the  identified  model  predict  these  measurements  better  than  a 
model  based  on  preflight  parameters? 

2.2  DATA  PROCESSING  ALGORITHMS 

A  large  number  of  methods  exist  for  performing  system  identification  data 
processing.  The  best  algorithm  for  any  given  application  depends  strongly  on 
the  type  of  model  and  on  the  nature  of  the  available  data.  No  one  type  of 
processing  algorithm  can  handle  all  possible  applications. 

This  section  outlines  three  processing  methods  (Table  2.1)  which  have 
been  found  to  be  effective  in  a  variety  of  applications.  These  methods  are: 

-  equation  error  minimization  methods, 

-  output  error  minimization  methods,  and 

-  simultaneous  state  and  parameter  estimation  methods. 

2.2.1  Equation  Error 

The  equation  error  minimization  methods  estimate  unknown  parameters  by 
choosing  them  to  minimize  a  performance  index.  A  continuous  dynamic  system 
must  be  represented  as: 


where  e  is  a  set  of  p  unknown  parameters  and  w  is  a  time-varying 
unobservable  disturbance.  An  analogous  formulation  exists  for  a  discrete 
dynamic  system.  The  performance  index  to  be  minimized  is: 

&  (&)  m  j  M  dx(t  ) 

e  -ftx(ti).  u(ti),ti,  e]  2  (1) 

The  equation  error  minimization  method  is  often  called  the  least  squares 
method  because  of  the  form  of  the  performance  index  ^e(e). 

The  effective  use  of  the  equation  error  minimization  requires  the  a 
priori  determination  of  system  states  x,  controls  u,  and  state 
derivatives  dx/dt  over  the  time  interval  of  the  test.  A  priori  here  means 
that  these  quantities  must  be  determined  before  the  unknown  system  parameters 
are  estimated.  This  determination  may  be  done  using  direct  measurements  or 
using  system  characteristics  which  are  independent  of  the  parameters.  For 
example,  an  unmeasured  state  derivative  may  be  determined  by  (very  carefully) 
numerically  differentiating  a  measured  state  time  history. 

The  term  w  is  a  stochastic  quantity  which  represents  unmeasurable 
process  disturbances  in  the  system,  w  includes  wind  gusts  and  unmodeled, 
high-order  aerodynamic  effects. 

The  special  advantage  of  the  equation  error  minimization  method  lies  in 
the  fact  that  many  nonlinear  dynamic  system  functions  f(x,  u,  t,  e)  are 
linear  in  the  parameters  e.  In  other  words, 

P 

f(x,  u,  t,  e)  -  z  e.f.(x,  u,  t)  +  fn+1(x,  u,  t)  (2) 

j-1  J  J  p  1 

The  functions  fj,  j*l,2,...p+l  are  independent  of  the  p  unknown 
parameters  ©j.  The  parameter  values  which  minimize  ^e(e)  can  be  found 
explicitly  using  linear  algebraic  operations  [14].  The  disadvantages  of  the 
equation  error  minimization  method  arise  primarily  from  the  requirement  for 
very  accurate  measurements  of  states  and  controls.  States  will  inevitably  be 
measured  with  some  error.  No  measurement  at  all  may  be  available  for  other 
states. 


System  Identification  Data  Processing  Methods 
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2.2.2  Output  Error 

Output  error  minimization  methods,  like  equation  error  minimization 
methods,  estimate  unknown  parameters  by  choosing  them  to  minimize  a 
performance  index.  The  dynamic  system  must  be  represented  as 

dx/dt  -  f(x,  u,  t,  e),  x(tQ)  »  g(e)  (3) 

y  »  h(x,  u,  t,  e)  +  v  (4) 

where  e  is  a  set  of  p  unknown  parameters  and  v  is  a  time-varying, 
unobservable,  additive  measurement  error.  The  performance  index  ^(e)  to 
be  minimized  is 

m  A  9 

2M  «  S  [yr  y(ti ,  e)r  (5) 

0  i*l  1  1 

Here  yi  is  the  observed  system  output  at  time  t..  y(t.,  e)  is  the 
system  output  y  predicted  for  time  t^  by  solving  the  system  state 
equations  and  measurement  equations  using  the  measured  system  inputs  u(t^) 
and  the  a  priori  parameter  values  e. 

The  effective  use  of  the  output  error  minimization  requires  the  very 
accurate  measurement  of  system  inputs  u  and  the  measurement  of  system 
outputs  y.  The  method  will  tolerate  errors  in  the  measurement  of  y. 

The  term  visa  stochastic  quantity  which  represents  instrument 
measurement  errors,  e.g.,  analog-to-digital  quantization  noise. 

The  special  advantage  of  the  output  error  minimization  method,  with 
respect  to  the  equation  error  method,  is  that  the  measurement  requirements  are 
greatly  relaxed.  The  method  does  not  require  the  accurate  measurement  of  all 
state  and  state  derivatives.  Rather,  it  is  effective  using  noisy  measurements 
of  the  limited  number  of  outputs  that  are  available. 

The  actual  determination  of  the  parameter  values  ©  which  minimize  the 
performance  index  j?Q(e)  is  computationally  more  complex  than  the 
minimization  3?e(©).  This  is  because  :?Q(e)  is  a  nonlinear  function  of 
the  parameter  set  ©.  Finding  the  minimizing  parameter  set  requires  an 
iterative  numerical  scheme  [15,  16].  The  application  of  such  numerical 
methods  is  often  not  straightforward. 
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The  principal  disadvantage  of  the  output  error  minimization  scheme  is 
that  it  does  not  explicitly  allow  for  the  presence  of  unmodeled  disturbances 
in  the  state  dynamics.  Such  disturbances  are  represented  in  the  equation 
error  method  by  the  term  w.  "Process  noise"  is  the  term  often  used  to 
describe  these  unmodeled  effects. 

It  should  be  noted  that  the  output  error  method  can  account  for  system 
dynamics  disturbances  of  unknown  magnitude  if  the  form  of  these  disturbances 
is  accurately  represented.  The  disturbance  w  must  be  explicitly  represented 
as 

w  *  w(t,  0)  (6) 

The  unknown  elements  of  the  disturbance  are  represented  using  part  of  the 
unknown  parameter  vector  0.  One  might  estimate  the  horizontal  plane 
components  of  a  steady  wind  present  during  a  flight  test,  for  example. 

2.2.3  Combined  State  and  Parameter  Estimation 

Methods  which  combine  state  and  parameter  estimation  are  required  if 
significant  levels  of  both  unknown,  unmodeled  disturbances  and  measurement 
errors  are  present  in  the  system  under  study.  The  performance  index  used  here 
is  very  similar  to  the  output  error  index  .2^(0).  However,  the  estimated 
outputs  y  are  now  direct  functions  of  the  observed  outputs  y.  The 
performance  index  is 

M  ? 

tfs(e)  -  2  [y.  -  y(t.,  0,  y)r  (7) 

s  i=l  1  1 

The  estimated  outputs  are  determined  using  both  the  system  dynamic  equations 
and  the  observed  values  of  the  outputs  themselves. 

Methods  for  the  determination  of  y  given  measurements  y  and  an 
assumed  form  of  the  system  dynamics  have  been  widely  studied  under  the  topics 
of  state  estimation  [17]  and  linear  system  observation  [18].  The  use  of  the 


Kalman  filter  to  estimate  y  in  the  modified  output  error  performance  index 
5?s  leads  to  "maximum  likelihood"  parameter  estimates  [19].  The  procedure 
requires  that  be  evaluated  as  a  function  of  ©  using  the  Kalman 

filter  to  estimate  y.  The  parameter  values  are  estimated  by  choosing  them  to 
minimize  ,2?s(e)  using  a  Gauss-Newton  method  [20].  The  use  of  this  maximum 
likelihood  estimation  procedure  often  allows  the  estimation  of  system  noise 
levels  as  well  as  of  parameters  describing  the  physical  plant. 

2.3  SOME  PRACTICAL  CONSIDERATIONS 

This  section  briefly  discusses  a  number  of  practical  considerations  which 
should  be  taken  into  account  in  an  effective  system  identification  data 
processing  scheme. 

2.3.1  Assumptions  Regarding  Measurement  Noise  Statistics 

A  common  problem  is  to  assume  that  the  measurement  error  should  be 
modeled  as  a  Gaussian  white  process  when  in  fact  systematic  errors  such  as 
bias  and  scale  factor  exist.  Systematic  measurement  errors  will  usually  cause 
larger  parameter  estimation  errors  than  random  noise  errors  of  the  same 
root-mean-square  level.  A  very  common  scale  factor  found  when  dealing  with 
any  instrument  using  electronic  pickoffs  is  -1.0.  This  is  due  to  simple 
polarity  errors  made  when  installing  the  instrument.  Reference  11  covers 
methods  of  assessing  the  relative  significance  of  systematic  measurement 
errors  and  random  measurement  errors. 


2.3.2  Number  of  Independent  Parameters  In  The  Model 

Problems  can  arise  from  an  attempt  to  fit  too  complex  a  model  to  the 

available  data.  The  chief  symptom  of  this  is  that  a  large  scatter  of 

estimated  parameter  values  will  be  seen  if  several  data  sets  are  used 

independently  to  estimate  values  for  the  same  parameter  set. 


2.3.3  Extrapolation  of  Results 


An  identified  model  should  not  be  used  to  predict  system  behavior  for 
operating  regimes  far  beyond  those  encountered  during  data  collection. 
Operating  regime  predictions  should  be  limited  in  both  input  bandwidth  and 
amplitude  to  those  tested. 

2.3.4  Excitation  of  All  System  Modes 

This  problem  can  be  avoided  by  careful  choice  of  inputs  during  the  test 
planning  stage.  A  second  solution  is  to  process  multiple  maneuvers 
simultaneously  which  contain  different  control  inputs.  By  doing  this,  the 
required  modal  information  is  extracted  from  a  set  of  simpler  maneuvers, 
rather  than  one  complicated  maneuver. 

2.3.5  Effective  Use  of  Sequential  Data  Processing  Schemes 

System  identification  data  processing  requires  the  computational  steps  of 
model  structure  determination,  parameter  estimation,  and  model  validation.  An 
additional  preliminary  step  of  prefiltering  measurements  may  also  be  required 
for  effective  use  of  an  equation  error  parameter  estimation  method.  An 
effective  overall  computational  scheme  may  require  that  the  operations  of 
prefiltering,  model  structure  determination,  and  parameter  estimation  be 
carried  out  in  a  sequential  rather  than  in  a  more  nearly  simultaneous  manner. 
Care  must  be  taken  to  ensure  that  the  algorithms  employed  at  any  given  stage 
do  not  remove  critical  information  from  the  data.  As  a  simple  example,  the 
bandpass  of  a  noise  prefilter  should  be  higher  than  that  of  the  modes  of  the 
system  to  be  identified. 

2.3.6  Process  Noise 


The  term  "process  noise"  refers  to  unmodeled  factors  in  the  state 
dynamics  of  the  system  being  identified.  Sources  of  process  noise  include: 
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(1)  unmeasured  environmental  disturbances  -  wind  gusts  acting  on  an 
aircraft,  for  example, 

(2)  unmodeled  nonlinearities  or  degrees  of  freedom  in  the  state 
dynamics,  and 

(3)  errors  in  measuring  input  signals. 

The  effect  of  process  noise  is  usually,  but  not  always,  to  degrade  estimation 
accuracies.  If  measurements  of  system  states  are  highly  accurate,  then  the 
process  noise  becomes  the  major  source  of  estimation  error.  Under  some 
circumstances,  process  noise  in  the  form  of  unmeasured  environmental 
disturbances  can  improve  estimation  accuracy.  The  environmental  disturbances 
might  excite  modes  of  the  system  which  are  not  excited  by  the  known  input  test 
signal . 

The  relative  significance  of  process  noise  in  an  identification  effort 
depends  roughly  upon  the  ratio 

r  «  RMS(process  noise) /RMS(known  inputs)  (8) 

where  RMS(  )  refers  to  the  root-mean-square  state  excursion  due  to  the 
indicated  source  of  excitation.  If  r  is  large,  then  the  process  noise  is 
significant.  If  r  is  small,  then  the  process  noise  is  not  significant.  It 
is  difficult,  however,  to  specify  a  value  of  r  indicating  the  boundary 
between  significant  or  insignificant  process  noise  levels  which  will  be  valid 
for  all  systems. 

Effective  system  identification  methods  exist  for  use  when  available  data 
contain  process  noise.  The  equation  error  formulation  is  preferred  if  all 
system  states  can  be  measured  or  estimated  accurately,  otherwise  the 
formulation  combining  state  and  parameter  estimation  will  be  required. 

2.3.7  Initialization 

Many  parameter  estimation  formulations  require  the  iterative,  numerical 
solution  of  nonlinear  equations.  The  output  error  and  the  combined  state  and 
parameter  estimation  formulations  fall  into  this  category.  Iterative 
numerical  algorithms  require  initial  estimates  of  parameter  values  in  order  to 


begin  execution  of  the  first  iteration.  Inaccurate  initial  estimates  may 
cause: 

(1)  convergence  of  the  estimation  method  (which  usually  employs 
some  form  of  performance  criterion  minimization  algorithm)  to  a 
local  minimum,  or 

(2)  divergence  of  the  estimated  parameter  values  as  iterations 
proceed.  Divergence  may  occur,  for  example,  if  the  values  of 
the  initial  parameter  estimates  cause  an  instability  in  the 
dynamic  system  model. 

An  effective  way  to  obtain  initial  parameter  estimates  for  starting 
iterative  algorithms  is  often  to  employ  the  equation  error  estimation 
formulation.  As  noted  in  Section  2.3.1,  the  equation  error  formulation 
usually  requires  only  the  solution  of  a  linear  set  of  algebraic  equations  in 
order  to  obtain  parameter  estimates.  Such  equations  may  be  solved  without  a 
priori  parameter  estimates.  An  effective  two-step  parameter  estimation 
procedure  is: 

(1)  estimate  initial  parameter  values  using  the  equation  error 
formulation,  then 

(2)  refine  these  estimates  using  either  the  output  error  or  the 
combined  state  and  parameter  estimation  formulations. 

The  values  of  parameters  estimated  using  the  equation  error  formulation  are 
sensitive  to  errors  in  measuring  states  (measurement  noise).  However,  the 
parameter  estimates  calculated  using  an  equation  error  criterion  even  with 
data  corrupted  by  measurement  noise  are  often  sufficiently  accurate  for  use  as 
start-up  values  for  iterative  algorithms. 

2.3.8  Numerical  Methods 


System  identification  algorithms  engender  a  variety  of  numerical 
mathematical  requirements.  Table  2.2  lists  four  of  these: 

-  solution  of  differential  equations, 

-  solution  of  linear  algebraic  systems  of  equations, 

-  solution  of  least  squares  problems,  and 

-  minimization  of  general  nonlinear  multivariable  functions. 
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Table  2.2 


Numerical  Methods  Used  in  System  Identification 


NUMERICAL  MATHEMATICAL  REQUIREMENT 

EFFECTIVE  fCTHOD  OF  APPROACH 

SOLUTION  OF  DIFFERENTIAL  EQUATIONS 

LINEAR 

TRANSITION  MATRIX  [23] 

NONLINEAR 

HULTISTEP  METHODS  (ADAMS-BASHFORTH)  [24] 

SOLUTION  OF  LINEAR  ALGEGRAIC 

SYSTEMS  OF  EQUATIONS 

POSITIVE  DEFINITE. 
SYMMETRIC 

CHOLESKY  [25] 

GENERAL 

GAUSSIAN  ELIMINATION  [25] 

SOLUTION  OF  LEAST  SQUARES  PROBLEMS 

LINEAR 

FACTORIZATION  OR  SQUARE  ROOT  METHODS  [26] 

NONLINEAR 

GAUSS-NEUTON  [27] 

MINIMIZATION  OF  GENERAL  NONLINEAR  MULTIVARIABLE  FUNCTIONS 

QUASI -NEWTON  [28] 

Effective  methods  to  handle  these  problems  range  from  the  classical  systematic 
elimination  method  of  Gauss  [25]  for  the  solution  of  systems  of  linear 
algebraic  equations  to  more  recent  developments  in  the  solution  of  linear 
least  squares  problems  [26]. 

An  important  consideration  in  using  any  of  the  methods  of  Table  2  is  that 
of  numerical  conditioning.  Numerical  conditioning  refers  to  the  sensitivity 
of  the  output  of  a  numerical  algorithm  to  small  changes  in  the  input  to  the 
algorithm.  For  example,  assume  that  in  solving  a  system  of  n  linear 
equations 

Ax  ■  b, 

the  matrix  A  is  known  exactly,  but  the  vector  b  is  subject  to 
uncertainty  6b.  The  norm  of  the  resulting  uncertainty  in  x, 

6x,  is  bounded  by  [25] 


where  uj  is  the  largest  eigenvalue  of  AAT  and  vn  is  the  smallest 
eigenvalue  of  AA^.  The  quantity 

cond(A)  =  Sviivn  (10) 

is  called  the  condition  number  of  A  and  is  always  greater  than  1.0.  Similar 
bounds  for  solution  sensitivities  exist  for  least  squares  parameter  estimation 
problems  [26]. 

The  condition  number  of  a  numerical  problem  can  give  insight  into  the 
precision  required  to  obtain  acceptable  accuracy  of  solution.  The 
uncertainty  «b,  for  example,  might  be  due  to  rounding  due  to  finite 
precision  in  computer  word  length.  If  the  condition  number  of  a  problem  is 
10®,  then  eight  significant  figures  of  accuracy  would  be  required  to 
maintain  1  percent  accuracy  in  the  solution. 
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2.4  CONSIDERATIONS  FOR  NONLINEAR  PARAMETER  IDENTIFICATION 


Methods  for  the  maximum  likelihood  identification  of  linear  state  dynamic 
models  are  well  established  [29].  Such  methods  have  been  applied  to  problems 
of  linear  modeling  of  aircraft  aerodynamics  using  flight  test  [22]  data. 

These  methods  are  sometimes  applied  in  a  piecewise  manner  to  fundamentally 
nonlinear  systems. 

The  intrinsic  nonlinear  nature  of  aircraft  aerodynamic  models  may  inhibit 
the  effective  use  of  linear  identification  methods.  For  example,  if  a  limited 
amount  of  data  is  available,  it  may  not  be  possible  to  identify  many  linear 
perturbation  models.  A  single  nonlinear  model  may  have  fewer  total  free 
parameters.  Also,  excursions  through  the  nonlinear  portion  of  the  model's 
dynamic  range  may  be  so  rapid  that  no  single  linearized  model  can  represent  a 
significant  portion  of  the  trajectory. 

The  use  of  a  nonlinear  model  may  be  required  if  the  goal  of  the  analysis 
is  to  determine  which  one  of  several  competitive  phenomenological  mathematical 
models  best  fits  available  data.  A  "phenomenological"  model  is  one  that  is 
constructed  from  fundamental  physical  principles.  Such  a  model  may  have  a 
very  complex  form  mathematically  but  may  have  a  minimum  number  of  unknown 
coefficients. 

There  are  certain  computational  difficulties  associated  with  the  use  of  a 
nonlinear  dynamic  model  in  a  maximum  likelihood  parameter  estimation  algorithm. 

Calculation  of  Sensitivities 

The  estimation  of  parameters  through  the  use  of  the  maximum  likelihood 
criterion  requires  the  maximization  of  the  likelihood  of  the  data  with  respect 
to  the  unknown  parameters.  The  determination  of  the  maximizing  parameter 
values  requires  numerical  optimization  techniques.  The  most  efficient  of 
these  [27]  are  descendents  of  the  Levenberg-Marquardt  nonlinear  least  squares 
method  [16,30].  These  algorithms  require  the  evaluation  of  the  partial 
derivative  of  modeling  residuals  with  respect  to  parameter  values.  This 
partial  derivative  is  often  called  a  "sensitivity".  The  calculation  of  these 
sensitivities  is  not  difficult  in  principle.  They  satisfy  differential 
equations  which  are  closely  related  to  the  system  dynamic  equations,  but  which 
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contain  terms  based  on  the  algebraic  partial  derivatives  of  the  dynamic 
equations  (see  Section  3.1).  The  difficulty  is  one  of  practice.  Any  time 
that  the  structure  of  the  nonlinear  model  is  changed,  then  the  sensitivity 
differential  equations  must  be  changed  also.  This  requires  tedious  algebraic 
differentiation  of  the  modified  dynamic  equations. 

Evaluation  of  Covariance  of  Parameter  Estimates 

A  parameter  covariance  matrix  can  be  estimated  using  the  Cramer-Rao  bound 
[31].  The  most  common  use  of  this  bound  assumes  that  the  errors  in  predicting 
the  response  of  the  system  are  due  to  an  additive,  white  (negligible 
autocorrelation)  random  process.  If  the  analyst  also  desires  confidence 
intervals  for  parameter  estimates,  then  the  additional  assumption  that  the 
errors  have  a  normal  distribution  must  also  be  made.  These  assumptions  are 
commonly  violated  when  a  nonlinear  system  is  modeled.  In  particular,  the 
whiteness  assumption  is  typically  violated. 

Generic  Model  Structure 

There  is  a  need  to  represent  nonlinear  functions  of  several  variables  in 
the  model  used  in  the  identification  algorithm.  Ideally,  a  single  generic 
form  should  represent  multidimensional  surfaces  of  arbitrary  shape.  These 
functions  represent  total  aerodynamic  force  or  moment  coefficients  as 
functions  of  angle  of  attack,  angle  of  sideslip,  and  angular  rates. 

2.4.1  Problem  Definition 


The  dynamic  system  is  modeled  as  nx  first  order  nonlinear  differential 
equations. 

1  =  f(x.»  u.  w,  t,  e)  (11) 

having  an  output  measured  at  discrete  times  tk 

y(tk)  =*  h  (x,  u,  tk,  &)  +  v  (tk)  (12) 
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Here 


x  *»  nx  component  state 

£  -  nu  component  inputs  measured  without  error 
&  -  nth  component  unknown  parameter  vector 

w  is  a  nw  component  random  input  (process  noise  source)  having 
statistics 


E[w  (ti)  wT(t^)]  »  Q( t i )  6^ 


The  scalar  v(tk)  *  vk  is  a  random  measurement  error  having  the 
statistics 

E(vk)  =  0  (15) 

E(v^)  *  r  (16) 

Note  that  the  assumption  of  scalar  measurements  does  not  cause  a  great  loss  of 
generality.  This  formulation  can  accomodate  multiple  sensors  simply  by 
assuming  that  the  interval  between  measurements  is  sometimes  very  small.  The 
only  loss  of  generality  regards  the  representation  of  correlation  of 
measurement  error  between  sensors. 

If  we  assume  that  the  stochastic  quantities  have  normal  distributions 
then  the  joint  density  or  likelihood  function  of  a  sample  of  nt 
measurements  y  (tk)  =  yk  is 


yj,  y2, 


»nt  :  2> 


exp  [-  \  (yk  -  yk(e))2/»k2(e)] 


An  extended  Kalman  filter  [17]  can  generate  both  the  measurement  estimates 
yk  and  the  measurement  uncertainties  oK.  The  estimate  of  e  having  the 


I 


smallest  variance  is  the  one  which  maximizes  ;  ©)  with  respect  to  e. 

The  maximization  of  SB  is  equivalent  to  the  minimization  of  the  negative  log 
likelihood  function  given  by 

-log  SB{^  ;  ©)  -  (18) 

\  l  j[yk  -  y  (»)]2/ok(£)  +  2  log  ok(e) 

K«l( 

Note  that  if  the  ak  are  known,  then  the  minimization  of  -log  SB  can  be 
treated  as  a  nonlinear  least  square  problem. 

Aspects  of  this  problem  addressed  here  are  the  following. 

(1)  The  minimization  of  Eq.  (18)  with  respect  to  ©  requires  the  use  of 
an  iterative  numerical  procedure  similar  to  a  Newton  or  quasi-Newton  method. 

If  the  ok  are  known,  then  the  most  effective  procedure  is  that  of 
Levenberg  and  Marquardt  [16,30].  A  drawback  to  the  Levenberg-Marquardt  method 
is  the  requirement  for  the  evaluation  of  ay/a©.  Direct  analog  finite 
difference  methods  [32]  avoid  this  problen  by  approximating  the  partial 
derivative  by  a  finite  difference.  Section  2.4.2  extends  the  finite 
difference  analog  of  the  Levenberg-Marquardt  procedure  to  the  more  general 
form  of  Eq.  (18). 

(2)  Validation  of  an  identified  model  should  include  the  determination 
of  confidence  intervals  or  variances  for  estimated  parameters.  These  may  be 
estimated  using  the  Cramer-Rao  bound,  which  states  that 


A  *  A  *  T  i 

£[(©-©)  (e  -  ©  )  ]  _>  M'1 

(19) 

where  M 

is  the  information  matrix,  given  by 

[M].^  .  -[32(log  SB{i  ;  ©)/a©1  3©^ 

(20) 

* 

and  © 

is  the  true  parameter  value. 

24 


Experienced  analysts  in  the  system  identification  field  are  aware  that 
the  Cramer-Rao  bound  is  usually  optimistic  [33],  That  is,  it  tends  to  predict 
parameter  variances  which  are  very  small  in  comparison  to  variances  observed 
among  estimates  derived  from  multiple  data  sets.  The  most  easily  implemented 
expressions  for  the  information  matrix  assume  that  the  measurement  errors  are 
not  autocorrelated.  If  autocorrelation  is  accounted  for,  then  the  Cramer-Rao 
bound  more  realistically  approximates  the  true  parameter  variance.  Section 

2.4.3  outlines  methods  based  on  principles  of  generalized  least  squares  which 
automatically  account  for  first  order  autocorrelation  of  measurements.  This 
method  produces  parameter  estimates  which  have  a  lower  variance  than  those 
which  do  not  account  for  autocorrelation. 

(3)  Modeling  nonlinear  aerodynamics  requires  the  representation  of 
nonlinear  functions  of  several  variables.  For  example,  pitch  moment  of  an 
aircraft  is  a  nonlinear  function  of  angle  of  attack  o,  angle  of  sideslip  s, 
pitch  rate  Q,  and  elevator  angle  «e. 

cm  '  cm<«  •  »•  «.  se>  <21> 

Much  work  in  the  identification  of  these  functions  has  been  based  on  their 
representation  as  multidimensional  polynomials  [34].  This  approach  is 
effective  for  local  models.  A  local  model  is  one  that  is  valid  over  a 
restricted  region  of  the  flight  envelope,  say  for  an  a  interval  of  10*. 

Such  models  often  use  expansions  for  functions  like  Cm  in  the  states  of 
degree  no  higher  than  two.  The  representation  of  a  global  model  using 
polynomial  expansions  may  require  very  high  order  polynomials.  The 
representation  of  a  pitch  moment  curve  for  one  particular  aircraft  through  a 
40*  angle  of  attack  region  requires  a  ninth  degree  polynomial  [35].  Section 

2.4.4  indicates  how  the  use  of  local,  low  degree  polynomial  models  leads  very 
naturally  to  a  spline  formulation  for  a  global  model. 


2.4.2  Derivative-Free  Minimization  of  the  Negative  Log  Likelihood  Function 

Finding  the  parameter  values  e  which  minimize  Eq.  (18)  requires  the  use 
of  a  numerical  optimization  scheme.  Typically,  each  iteration  of  the 
algorithm  updates  e  as 


Vl  -  2n  + 

where  ae  satisfies 


Mae  » 


g.  =  alog  SB  (^,  e)/ae. 

and  M  is  defined  in  Eq.  (20). 

Differentiating  Eq.  (18)  gives  [35] 


9i  *  *  k£1{Cvk(avk/a®i)  ‘  7  \(aok/3ei)/ok 


4  °k  /39i]/ok} 


/  2  2 
-  £  j  (3'»k/3©i)(8'»lc/aej)/ok+  [-(av^/ae^tao^/aej) 


-(3vk/a©j)(a ak/ae.)  +  vk(aek  /a®i )(ae^/aej.) 


+  °k^3^'Jk  /3®i39j)  *  vk(a2ok/aeiaeJ.)]  vk/oJ 


*j(H2  ,k2/39)  S9j)/,2} 
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where 


vk  "  *k  -  *k  <£>• 


(27) 


Following  the  spirit  of  the  Levenberg-Marquardt  method,  we  simplify  the 

expression  for  the  Hessian  by  dropping  the  terms  proportional  to  If 

the  model  fits  the  data  well,  vk  should  approach  zero  near  convergence. 

With  the  exception  of  the  last  term  of  Eq.  (26),  it  is  now  possible  to 

evaluate  both  gradient  and  Hessian  if  only  the  first-order  sensitivities 
2 

and  ok  to  changes  in  e  are  known. 

It  is  possible  to  derive  analytically  ordinary  differential  equations 

2 

for  avk/3©i  and  ao^/ae^  When  the  plant  dynamic  equations  are 
linear,  these  differential  sensitivity  equations  (also  called  “sensitivity 
equations")  have  a  particularly  simple  form  [36].  When  the  plant  dynamics  are 
nonlinear,  however,  a  much  more  practical  method  is  to  approximate  the  partial 
derivatives  with  finite  differences. 

aVfc/aSi  *  tv|((®  +.£j)  ~  vl<(±)]/!.£.il  (28) 

or 

avfc/aoi  =  [vk(e  +  £})  -  vk(e  -  ei)]/2|£il  (29) 

where  U-j)j  =  c6ij* 

2 

Partial  derivatives  of  o  (e)  are  approximated  similarly. 

The  "direct  analog"  [32]  type  of  optimization  algorithms  use  the  finite 
difference  approximations  to  the  partial  derivatives  as  direct  substitutes  for 
the  partial  derivatives.  These  algorithms  have  been  studied  carefully  for  the 
solution  of  the  nonlinear  least  square  problem  and  have  been  found  to  have 
convergence  properties  nearly  identical  to  algorithms  which  use  analytic 
expressions  for  the  sensitivities. 


Alternate  methods  for  minimization  of  the  negative  log  likelihood 
function  without  evaluating  derivatives  use  algorithms  which  are  not  direct 
analogs  of  derivative  methods.  We  have  tested  one  of  these  methods  [37]  in  a 
nonlinear  system  identification  algorithm  and  found  it  to  be  not  as  effective 
as  the  direct  analog  method. 

In  partial  implementations  of  the  finite  difference  method  for  the 
solution  of  system  identification  problems,  we  typically  simultaneously  solve 
a  set  of  nx  nominal  state  equations  together  with  nth  sets  of  nx  perturbed 
parameter  state  equations.  This  is  required  in  order  to  implement  the 
one-sided  approximation  to  the  sensitivity  given  in  Eq.  (28). 

The  last  term  of  Eq.  (26)  cannot  be  eliminated  by  assuming  that  vk  is 

small  near  the  minimum.  The  term  cannot  be  constructed  from  first  order 

2 

partial  derivatives  of  v  and  a  .  It  can  be  estimated,  however,  using 
methods  similar  to  those  employed  in  solving  large  residual  nonlinear  least 
square  problems  [38].  This  topic  will  not  be  treated  in  further  detail  here. 

2.4.3  Autocorrelated  Measurement  Errors 

It  is  desireable  to  determine  the  expected  accuracy  of  the  parameters 
estimated  from  flight  data.  The  usual  method  for  doing  this  is  to  compute  a 
parameter  covariance  matrix  as  the  inverse  of  the  Fisher  information  matrix 
M. 

Methods  of  generalized  least  squares  [39]  indicate  expressions  for 
parameter  estimation  covariance  using  the  assumption  that  measurement  errors 
are  autocorrelated.  This  autocorrelated  process  is  the  output  of  a  first 
order  difference  equation  driven  by  white  noise.  Not  only  parameter 
covariances  but  also  parameter  estimates  themselves  are  altered  by  the 
autocorrelation  assumption.  An  estimation  algorithm  which  does  not  explicitly 
account  for  the  measurement  error  autocorrelation  will  still  produce  unbiased 
parameter  estimates.  However  the  actual  variance  of  such  estimates,  as 
opposed  to  the  Cramer-Rao  predicted  variance,  will  be  higher  than  those 
produced  by  an  algorithm  which  does  explicitly  account  for  autocorrelation. 


I  E 

i  p 

B 


i  r 

3  E 

«  L 

C  - 


We  consider  here  only  the  output  error  case  (no  process  noise)  of  the 
maximum  likelihood  estimator  for  dynamic  systems.  The  estimated  measurements 
y  are  functions  of  the  unknown  parameter  set  e.  The  information  matrix  M 
is  given  by 

T  (30) 

a 

where  J  is  the  matrix  of  sensitivities 

[J]k2  "  ”  (31) 

for  the  case  of  white  measurement  noise.  Iterations  of  the  identification 
algorithm  solve 

Mae  =  -£  (32) 

where 

I  =  7  JT  «  (33) 

a 

k  -y(tk)  -y(tk> 

The  covariance  of  the  measurement  error  is  a  diagonal  matrix 

E(vvT)«o2I  (34) 

where  vk  »  v(tk)  *  the  measurement  error  at  tfc.  (35) 

Now  suppose  that  the  measurement  errors  have  a  nondiagonal  covariance 
matrix  of 

E[vvT]=o2V  (36) 

For  a  purely  linear  estimation  problem,  (i.e.  y  *  Je  +  o),  the  "general ized 
least  squares"  estimate  of  e  satisfies  [39] 

[JTV-1J]e  -  jV1*  (37) 

The  covariance  of  the  parameter  estimates  is 

E[(£  -  e*)(e  -  e*)T]  =  a2(jV1J)'1  =  Var(e)  (38) 


c 


* 
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Any  other  linear,  unbiased  estimator  has  a  covariance  matrix  which  exceeds 
that  given  in  Eq.  (38). 

If  the  noise  vector  v^  is  generated  by  a  first  order  autoregressive 
process 


v(tk)  =  pvU^)  +  ek  (39) 

where  ek  is  a  zero  mean,  constant  variance  process,  then  V  has  the  form 


The  generalized  least  square  estimator  for  ©  can  be  easily  implemented  by 
writing  Eq.  (37)  as 

C(P  J)t(p  J)]'1  i=(P  J)tpi 

The  multiplications  P  J  and  P  ^  are  simple  because  £  is  sparce. 


(42) 


i 


a  r 


For  application  to  the  output  error  system  identification  problem,  a 
nonlinear  least  square  problem,  Eq.  (42)  is  applied  to  parameter  variations 
ae  on  each  iteration  of  a  successive  approximation  algorithm. 


(P  J)T(P  J)ae  =  (P  J)  T  P  -  £  (»)] 


If  p  is  unknown,  then  it  can  be  estimated  using 


~  nt  I 

E  vZ(tJ 
k-1  K 


nt  -  1 
nt  -  nth 


Each  iteration  calculates  a  Ae  value  using  P  evaluated  by  Eq.  (44),  with 

v(tk)  calculated  from  e  at  the  end  of  the  previous  iteration. 

o 

If  the  measurement  noise,  o  ,  is  unknown,  it  can  be  estimated  using 


32  =  Tt'“nth(P',)T(Pv)  * 


The  covariance  of  the  parameter  estimates  is 


Var(£)  ■  TiTTnth  (p^T(Pv)[(PJ)T(PJ)3~1 


where 


P  =  9 A  -  p2  *  ^ 

2.4.4  Spline  Model  Structure 

The  determination  of  a  nonlinear,  quasi-static  aerodynamic  (or 
hydrodynamic)  model  requires  definition  of  a  coefficient  function  having  a 
general  form 

C  =  C(a,  e,  u',  Rg,  Fr>  M)  (* 
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where  a  and  s  are  relative  flow  angles,  w'  is  a  dimensionless  angular 
rate  vector,  is  a  control  vector,  and  Rg,  F  ,  and  M  are  the 
dimensionless  numbers  of  Reynolds,  Froude,  and  Mach.  The  model  structure 
determination  problem  for  identification  of  aerodynamic  models  usually  refers 
to  the  problem  of  determining  a  mathematical  form  for  this  multivariable 
function. 

Spline  functions  are  effective  ways  to  represent  these  coefficient 
functions.  A  one  dimensional  spline  function  is  a  piecewise  polynomial 
function  having  certain  continuity  conditions  between  pieces.  Figure  2.2 
illustrates  a  one  dimensional  cubic  spline.  C(a)  might  represent  pitch 
moment  as  a  function  of  angle  of  attack.  C(o)  here  is  a  cubic  polynomial  on 
each  of  the  three  regions  indicated.  The  function  is  everywhere  continuous 

and  has  continuous  first  and  second  derivatives.  The  points  aj,  a^, 

a^,  a4  are  called  the  knots  of  the  spline. 

The  spline  function  has  several  properties  which  make  it  an  effective 

interpolating  function  [40]. 

(1)  The  spline  in  Figure  2.2,  for  example,  is  uniquely  determined  once 
the  values  of  the  function  at  the  four  knots  are  known  and  certain  end 
conditions  are  specified. 

(2)  The  shape  of  the  interpolating  function  is  not  overly  sensitive  to 
the  function  values  at  the  knots.  Small  changes  in  these  values  do  not  cause 
overly  large  changes  in  interpolated  function  values  between  knots. 

(3)  The  interpolating  function  C(a)  has  an  optimal  smoothness 
property.  It  is  the  unique  function  which  interpolates  the  specific  values  at 
the  knots,  has  the  continuity  conditions  listed  above,  and  has  the  minimum 
mean  square  curvature. 


Spline  function  representation  of  nonlinear  aerodynamic  or  hydrodynamic 
coefficient  functions  may  be  readily  identified  using  either  maximum 
likelihood,  equation  error  or  output  error  techniques.  The  analyst  must 
specify  the  number  of  knots  and  their  locations.  The  parameters  to  be 
identified  are  then  the  function  values  at  the  knots.  The  identification  of 
the  C(a)  curve  in  Figure  2.2  would  require  the  estimation  of  four  parameters. 

In  Figure  2.2,  the  coefficients  -  C4  are  the  function  values  at 
the  knot  locations  -  a^.  These  coefficients  are  the  parameters  which 
will  be  estimated  by  the  identification  algorithm.  The  piecewise  cubic 
polynomials  provide  cubic  interpolation  of  for 

a  in  the  range  [a^,  a^]  and  linear  extrapolation  for  a  outside  this 
range. 


of 


Each  of  the  piecewise  cubic  polynomials  is  defined  over  the  entire  range 
[-ox-  -  <  +  «].  The  function  C(a)  is  a  linear  combination  of  the  K, 


Ai 


basis  functions. 


'At 


are  such  that 


This  makes  the  coefficients  in  the  linear  combination  equal  to  the  Ci 


values.  The  array  in  Figure  2.2  defines  the  Kft.  functions  over  the  five 


a  regions. 


The  identification  of  spline  functions  is  most  effective  when  used  with 
derivative  free  methods  to  minimize  the  negative  log  likelihood  function. 

Such  methods  do  not  require  the  explicit  calculation  of  the  sensitivity  of  the 
spline  function  to  changes  in  the  parameters  which  define  the  spline.  The 
only  requirement  is  for  the  evaluation  of  the  spline  coefficients  given 
function  values  at  the  knots,  and  for  the  evaluation  of  the  function  at 
intermediate  points  given  the  spline  coefficients.  Each  iteration  of  the 
“direct  analog"  method  requires  the  evaluation  of  the  innovations  (tk)  for 
nominal  and  for  perturbed  parameter  values. 


Methods  exist  for  the  use  of  multidimensional  spline  functions  to 
represent  smooth  surfaces  [41].  Intermediate  methods  are  also  useful.  An 
intermediate  method  represents  variation  of  a  function  in  one  dimension  with  a 
spline  function  and  representation  in  other  dimensions  with  other  types  of 
functions,  such  as  low  order  polynomials. 
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Figure  2.2  Cubic  Interpolating  Polynomials 
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2.5  MODELING  APPROACH 


The  approach  selected  for  modeling  nonlinear  aerodynamic  characteristics 
produces  system  identification  results  that  can  be  used  to  validate  preflight 
estimated  aerodynamic  models.  The  models  are  used  for  flight  simulators  and 
for  making  predictions  of  aircraft  performance,  stability  and  control 
characteristics.  These  aerodynamic  models  must  account  for  the  effect  of  a 
number  of  flight  condition  and  aircraft  configuration  variables.  The  ''art"  in 
formulating  the  models  is  to  represent  the  total  aerodynamic  coefficient  by  an 
incremental  buildup,  with  each  increment  described  by  one  or  two  independent 
variables.  This  process  is  illustrated  by  the  following  example  for  the 
rolling  moment  coefficient  equation. 

(1)  Select  the  independent  variables: 

=  f(a,  6,  P,  R,  4^,  4^) 

i.e.,  rolling  moment  coefficient  is  a  function  of  angle  of  attack, 
sideslip,  roll  and  yaw  rate,  rudder  and  aileron  position. 

(2)  Partition  independent  variables  into  reasonable  groups: 

C  *  aC  +  aC  +  aC 

*  “'SIDESLIP  “’DYNAMIC  “RUDDER 


AILERON 

(3)  Select  functional  relationships  for  each  group: 

aC  =  f(a,  8) 

“SIDESLIP 

■  C.  (a)  •  6 
8 


DYNAMIC 


f(a,  P)  ♦  f (a,  R) 


C*  (a)(Pbw/2VT)  +  Cz  <a)(Rbw/2VT) 
P  R 


» 
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iCp  *  ^p) 

RUDOER  K 


Cj  (<*)  *  5r 


*  f(«,  fi/v) 

»  ti  rnA»i  n 


AILERON 


Co  (<*)  *  6. 


For  this  model  formulation,  each  of  the  stability  derivatives,  (i.e., 

C  )  is  modeled  as  a  nonlinear  function  of  angle  of  attack.  By  using  a 

cubic  interpolation  spline,  as  described  in  Section  2.4.4,  the  parameter 

identification  algorithm  solves  for  C  at  specific  values  of  angle  of 

attack  (i.e.,  the  knots  of  the  spline).  This  procedure  is  illustrated  in 

Figure  2.3  which  shows  the  identified  variation  of  C  with  angle  of  attack. 

For  the  F-4S  parameter  identification  study,  was  identified  for  a  = 

6 

5*.  15*,  25*.  and  35*.  The  lower  four  parts  of  Figure  2.3  illustrate  the 
variation  of  the  interpolation  splines  with  angle  of  attack. 


Because  these  interpolation  splines  are  scaled  by  the  appropriate  value 

of  Cz  (i.e.,  the  a  -  5*  spline  is  scaled  by  the  value  of  Cz  for 

.  e 

o  ■  5*),  the  summation  of  the  four  interpolation  splines  defines  the  value  of 

C  for  any  value  of  angle  of  attack.  It  should  also  be  noted  that  each 
b 

interpolation  spline  has  the  value  of  C.  when  a  equals  its  knot 

0 

value,  and  it  is  zero  for  other  knot  values. 

The  spline  formulation  is  suitable  also  for  representing  installed 
propulsion  system  performance  models  and  test  instrumentation  calibration 
factors.  More  complicated  models  can  be  represented  by  a  bicubic  spline 
formulation. 
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III.  FLIGHT  TEST  PROGRAM  DESCRIPTION 


This  section  describes  the  test  aircraft,  the  flight  test  instrumentation 
system,  and  the  test  conditions  evaluated  for  the  F-4J  identification  program 
described  by  this  report.  The  descriptions  presented  in  this  section  are 
provided  as  background  information  for  the  discussions  on  model  structure 
determination  and  parameter  identification.  This  material  is  from  Refs.  42 
and  43. 

3.1  TEST  AIRCRAFT  DESCRIPTION 

The  test  aircraft,  F-4S  Bu  No.  157286  (McAir  No.  3846)  was  modified  under 
McAir  ECP  1049  Rl,  "F-4S  Two  Position  Leading  Edge  Slat  Retrofit  Program." 
Major  changes  to  distinguish  this  slatted  F-4S  model  from  previous  F-4S 
aircraft  are  noted  as  follows: 

(1)  two  position  leading-edge  slats; 

(2)  movable  inboard  leading-edge  flaps; 

(3)  42-unit  angle-of -attack  system; 

(4)  incorporated  12-pound  overbalance  weight  to  longitudinal 
control  system;  and 

(5)  J79-GE-10B  smokeless  combustor  engines. 

The  basic  motivation  for  adding  leading-edge  devices  to  the  F-4  was  to 
improve  its  high  angle  of  attack  aerodynamic  characteristics.  By  forestalling 
separated  wing  flow  with  the  leading-edge  devices,  it  would  be  possible  to 
expand  the  aircraft's  maneuvering  envelope  and  improve  its  handling 
qualities.  Figure  3.1  illustrates  some  of  the  aerodynamic  benefits  of 
leading-edge  devices  by  comparing  slat  on/off  wind  tunnel  data  [42]. 

The  following  subsections  present  a  general  description  of  the  aircraft 
configuration,  a  tabulation  of  pertinent  geometric  characteristics,  and  mass 
properties. 
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3.1.1  Aircraft  Configuration 


The  slat  equipped  F-4J  is  a  two-place,  all-weather  fighter  capable  of 
performing  as  a  missle-launching  or  cannon-firing  interceptor  or  an 
intermediate  and  long-range  attack  bomber.  General  configurational  features 
of  the  F-4S  are  illustrated  by  an  in-flight  photograph  (Figure  3.2)  of  the 
test  aircraft.  The  aircraft  is  powered  by  two  General  Electric  J79-10B 
turbojet  engines  with  automatically  controlled  external  compression  air 
inlets.  The  basic  design  is  characterized  by  a  low  aspect  ratio  wing  swept 
back  45*  at  the  25  percent  chord  line,  and  an  all-movable  slotted  leading-edge 
stabilator  with  23-1/4*  of  negative  dihedral.  Lateral  control  is  provided  by 
a  spoiler-aileron  combination.  Directional  stability  and  control  are  provided 
by  a  vertical  fin-rudder  combination.  The  following  paragraphs  provide 
functional  description  of  the  wing  leading  and  trailing-edge  high  lift 
devices. 

The  F-4S  flap/ si  at  system  is  an  integrated  system  that  provides  for 
automatic  slat  configuration  changes  during  in-flight  maneuvering  and 
selective  flap/slat  configurations  for  take-off  and  landing.  Each  wing  has 
two  leading-edge  slats  (one  on  the  inboard  wing  panel  and  one  on  the  outboard 
folding-wing  panel),  one  leading-edge  flap  mounted  adjacent  to  the  fuselage, 
one  trailing-edge  flap  which  incorporates  boundary  layer  control,  and  one 
aileron  with  droop  capability  (Figure  3.3).  The  entire  system  is  electrically 
selected  through  solenoid-operated  selector  valves  and  hydraulically  actuated 
using  the  utility  hydraulic  system. 

The  leading-edge  slats  are  designed  to  operate  a  two-position  system: 
slats  out  (extended)  for  the  high-lift,  maneuvering  configuration  and  slats  in 
(retracted)  for  the  clean-cruise  configuration.  The  slats  consist  of  inner 
and  outer  sections,  each  outer  section  having  two  segments  to  minimize  the 
effects  of  outer  wing  bending.  Both  segments  of  the  outer  section  are 
positioned  by  a  single  hydraulic  actuator  and  the  inner  section  is  positioned 
by  another  hydraulic  actuator.  All  four  slat  actuators  are  powered  by  the 
utility  hydraulic  system  and  all  slat  sections  move  simultaneously  when 
actuated.  The  inner  slat  section  translates  between  the  retracted  and 


extended  positions  and  the  outer  slat  section  pivots  between  the  retracted  and 
extended  positions.  All  slat  sections  are  mechanically  locked  in  both 
retracted  and  extended  positions  by  overcenter  linkages. 

The  leading-edge  flaps  are  designed  to  operate  as  a  two-position  system: 
flaps  down  for  take-off,  slow  flight,  and  landing;  and  flaps  up  for  the  clean- 
cruise  configuration.  Each  leading-edge  flap  is  positioned  by  a  hydraulic 
actuator.  Each  flap  is  held  in  the  extended  position  by  hydraulic  pressure 
and  locked  in  the  retracted  position  by  overcenter  linkages.  Leading-edge 
boundary  layer  control  has  been  removed  from  the  slat-modified  F-4J. 

The  trai ling-edge  flaps  are  designed  to  operate  as  a  three-position 
system:  flaps  up  for  the  clean-cruise  configuration;  flaps  one-half  for  the 
take-off  configuration;  and  flaps  down  for  the  slow  flight  and  landing 
configuratiion.  Each  trailing-edge  flap  is  positioned  by  a  hydraulic 
actuator.  Each  flap  is  held  in  the  extended  position  by  hydraulic  pressure 
and  locked  in  the  retracted  position  by  internal  locks  in  the  hydraulic 
actuating  cyliner.  Trailing-edge  flap  boundary  layer  control  is  operative 
only  when  the  flaps  are  in  the  full  down  position. 

Each  aileron  droops  (deflects  16-1/2*  down)  whenever  one-half  or  full 
flaps  are  selected.  This  is  accomplished  by  utilizing  an  aileron  droop 
cylinder  which  is  repositioned  by  an  electromechanical  droop  aileron 
actuator.  Both  aileron  systems  continue  to  function  as  originally  designed 
when  drooped,  except  that  the  aileron  neutral  point  is  16-1/2*  lower. 

Flap/slat  control  is  provided  by  the  existing  switch  in  the  front 
cockpit.  The  switch  is  airfoil  in  shape  and  has  newly  labeled  positions,  from 
top  to  bottom,  of  UP-NORM,  1/2-OUT,  and  DN-OUT.  When  the  switch  is  in  the 
UP-NORM  position,  all  flaps  are  fully  raised  and  the  slats  operate 
automatically  as  a  function  of  angle  of  attack  (provided  that  aircraft  weight 
Is  off  the  landing  gear).  As  long  as  aircraft  weight  is  off  the  gear,  the 
slats  are  extended  when  AOA  exceeds  11.5  units  and  retracted  when  AOA 
decreases  to  10.5  units.  With  aircraft  weight  on  the  landing  gear,  the  slats 
retract  if  UP-NORM  is  selected  regardless  of  AOA  value  as  a  safety  feature  (to 
prohibit  inadvertent  extension  or  retraction  of  the  slats  on  the  deck  due  to 
AOA  vane  movement). 


Two  airspeed  retraction  switches  are  provided  to  protect  the  slats, 
flaps,  and  ailerons  from  structural  damage  if  they  are  inadvertently  left 
extended  above  their  structural  airspeed  limits.  The  slats  airspeed  switch 
retracts  the  slats  when  accelerating  through  an  airspeed  of  approximately  585 
KCAS.  The  flaps'  airspeed  switch  retracts  the  flaps  and  ailerons  when 
accelerating  through  an  airspeed  of  approximately  237  KCAS.  If  the  slats  or 
flaps  are  retracted  by  their  respective  airspeed  switches,  they  will  return  to 
the  selected  flap/slat  switch  position  when  the  airspeed  is  reduced  below 
specific  limits. 

3.1.2  Aircraft  Geometry 

General  information  on  the  geometric  and  mass  characteristics  of  the  test 
aircraft  are  presented  in  this  section.  These  data  are  required  for  the 
solution  of  the  equations  of  motion  and  in  the  reconstruction  of 
nondimensional  aerodynamic  force  and  moment  coefficients  from  flight 
measurements.  Dimensions  of  the  F-4S  are  shown  in  Figure  3.3,  and  the 
required  reference  geometry  is  defined  below. 

-  wing  area,  sw  *.  530  ft2 

-  wing  span,  bw  .  38.67  ft 

-  wing  mean  aerodynamic  chord,  c"w  =  16.04  ft 

Moment  arms  for  calculating  thrust  moments  and  the  displacements  between 
the  c.g.  and  sensor  locations  are  presented.  The  latter  set  of  geometry  is 
required  to  reconstruct  equivalent  c.g.  measurements  since  the  equations  of 
motion  are  defined  in  terms  of  c.g.  quantities.  (The  equations  which  are  used 
to  define  equivalent  c.g.  measurements  are  presented  in  Section  IV.)  The 
location  of  the  c.g.,  sensors,  and  the  thrust  reaction  point  are  expressed  in 
terms  of  a  fuselage  station  (FS),  a  water  line  (WL),  and  a  butt  line  (BL)  in 
inches.  Moment  arms  are  then  defined  in  terms  of  axial  (£x),  lateral 
Uy),  and  vertical  displacements  (£z)  from  the  c.g.  in  feet.  The 
following  sign  convention  is  used:  tx  >  0  ahead  of  the  c.g.,  >  to  the 

right  of  the  c.g.,  and  £z  >  0  below  the  c.g.  The  following  set  of  data 
defines  the  c.g.  location  and  the  location  and  moment  arms  for  the 


accelerometer  package,  nose  boom  a  and  0  vanes  and  the  net  thrust  reaction 
point. 


c.g.  Location 

FS  »  261  +  1.92  XCg  (where  XCg  ■  c.g.  location  in  percent  c^) 
WL  -  28 

BL  -  0 


Accelerometer  Package 
FS  -  347.62 
WL  -  12.75 
BL  *  4.25  (right) 


-  7.22  +  .16  X 

eg 


.35 


l  -  1.27 
ZACC 


Noseboom  a  Vane 


FS  - 
WL  - 
BL  - 


-  65.5 

8.02 

-  8.5  (left) 

-  27.21  +  16  X„„ 
eg 


-  -  .71 


lz  «  1.67 
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Noseboom  g  Vane 
FS  =  -  70.15 


WL  =  -  .48 


BL  =  0 


27.60  +  16  X, 


=  -  2.38 


Net  Thrust  Reaction  Point  (Used  for  Phase  I) 
FS  -  313.5 
WL  =  32.4 

BL  =  23.8  left  and  right 
l  -  -  4.38  +  .16  X_ 


±  1.98 


-  .37 


Gross  Thrust  Reaction  Point  (Phases  II  and  III) 
FS  »  518 
WL  *  13 
BL  «  *23 

Inlet  Ram  Force  Reaction  Point 
FS  »  162 
WL  -  30 


BL  -  *35 


3.2  TEST  INSTRUMENTATION 


The  airborne  data  acquisition  system  used  for  the  flight  test  program 
included  a  three-axis  rate  gyro,  a  vertical  gyro,  a  directional  gyro,  engine 
RPM  and  fuel  flow  measurements,  a  three-axis  linear  accelerometer,  position 
and  force  transducers,  and  a  test  airdata  system.  The  airdata  system  uses  a 
noseboom  which  has  a  pitot-static  head  for  the  measurement  of  impact  pressur, 
static  pressure  and  temperature.  The  noseboom  also  has  vanes  for  measurements 
of  angle  of  attack  and  sideslip. 

Figure  3.4  shows  test  aircraft  body  axis  inertias  and  product  of  inertia 
as  a  function  of  gross  weight.  These  data  have  been  derived  for 

-  zero  fuel  weight  *  33,243  lbs 

-  no  armament 

-  tank  7  fuel 

-  flight  test  pod 

-  gear  up 

3.3  TEST  FLIGHT  CONDITIONS 
3.3.1  Phase  I 


Ten  of  the  Phase  I  flight  conditions  were  selected  for  system 
identification  processing.  These  flight  conditions  are  listed  in  Table  3.1 
and  their  selection  was  based  on  having  one  test  condition  for  each  type  of 
control  input.  Aircraft  weight  and  center  of  gravity  location,  run  time,  and 
initial-minimum-maximum  values  for  pressure  altitude,  equivalent  airspeed, 

Mach  number,  and  angle  of  attack  are  tabulated  in  Table  3.2  for  each  flight 
condition.  All  test  conditions  are  for  subsonic  flight  (M  <  .6)  along  the 
minimum  speed  boundary  of  the  operational  envelope. 

When  the  Phase  I  test  conditions  were  flown,  a  number  of  test  points  were 
flown  back  to  back.  After  each  run,  the  pilot  dropped  the  nose  to  regain 
airspeed  and  then  retrimmed  for  the  next  test  point.  Because  the  nose  dropped 


Table  3.1 

Flight  Conditions  Processed  for  Model  Structure  Determination 


FLIGHT  TEST 
NUMBER 

TEST  INPUT 

COMMENTS 

1004-5 

Full  Aft  Stick 

“MAX  47*1  Moderate 
Pitch-up  Rate 

1004-11 

Full  Aft  Stick  with  Simultaneous 
Lateral  Stick  Doublet 

“MAX  52';  Rapid 

Pitch-up  Rate;  Sustained 
Oscillation  for  a  >  40* 

1005-7 

Moderate  Aft  Stick  with 

Continued  Series  of  Lateral 

Stick  Doublets 

Oscillatory  a  in  15*  <  a 
<  40*  Range 

1005-10 

Moderate  Aft  Stick  with 
Simultaneous  Pedal  Doublet 

Oscillatory  a  In  15*  <  a 
<  30*  Range 

1006-4 

Full  Aft  Stick  at  Slow  Applica¬ 
tion  Rate  with  1-1/2  Pedal 

Doublets 

“MAX  ~  53* 

1006-5 

Moderate  Stick  with  Repeated 

Pedal  Doublets 

Oscillatory  a  in  15*  <  a 
<  34* 

1006-8 

Simultaneous  3-Axis  Doublets 
with  Repeated  Lateral-Directional 
Doublets  at  end  of  Maneuver 

Oscillatory  a  In  9*  <  a 
<  34*  Range 

1007-5 

Full  Aft  Stick  with  Full 

Lateral  Stick  2  Sec  Command 

“MAX  ~  42*;  Moderate 
Pitch-up  Rate 

1009-5 

Full  Aft  Stick  with  4  Sec 

Lateral  Command 

“MAX  ~  57*;  Moderate 
Pitch-up  Rate 

1010-7 

Full  Aft  Stick  with  2  Sec 

Pedal  Command 

“MAX  W55’;  Moderate 
Pitch-up  Rate 

[•] 


Condition  Summary  for  Phase  I  Flight  Test 


through  vertical  and  the  airplane  rolled  over  on  its  back  during  the  recovery 
maneuver,  this  phase  of  the  flight  record  was  not  used  for  system 
identification.  When  data  files  were  established  for  each  flight,  the  final 
time  was  selected  by  one  of  the  following  conditions: 

-  |  0 1  <  120* 

-  e  <  -20* 

-  stick  position  has  been  pulled  beyond  «s  >  10*  and  then  pushed 
through  trim 

3.3.2  Phase  II 

Data  for  the  second  phase  were  obtained  from  a  flight  which  took  place  on 
27  March  1980.  The  flight  test  plan  for  this  flight  (W/0-7)  was  prepared 
jointly  by  NATO  and  SCT  personnel. 

The  flight  test  program  was  conducted  in  a  clean  configuration  with  the 
throttles  fixed  for  each  maneuver  and  with  flight  near  and  beyond  stall.  Test 
conditions  were  generally  initiated  from  wing  level,  constant  altitude  flight 
with  a  -  10*  and  M  %  .6.  For  some  test  conditions,  the  pilot-applied  variable 
aft  stick  to  maneuver  the  aircraft  through  the  desired  angle-of -attack  test 
range.  For  some  of  the  stall  entry  conditions,  lateral  stick  and/or  pedal 
doublets  were  combined  with  the  longitudinal  stick  command.  Other  flight 
conditions  included  single-axis  and  multi-axis  sequenced  doublet  inputs. 

These  inputs  were  made  by  the  pilot  and  were  not  intended  to  be  repeatable  nor 
precise  with  regard  to  their  spectral  characteristics,  but  were  generally 
effective  in  creating  large-amplitude  motions.  The  overall  test  goal  was  to 
force  the  aircraft  through  a  broad  range  of  test  conditions. 

A  summary  of  flight  test  conditions  from  flight  W/D-7,  which  were  used 
for  the  identification  study,  is  presented  in  Table  3.3. 

3.3.3  Phase  III 

Three  flights  in  the  F-4S  were  flown  for  Phase  III  in  order  to  update  the 
existing  F-4S  simulation  data  base.  Test  techniques  used  include  variable 
frequency  sweep  inputs  in  all  three  axes  for  small  amplitude  maneuver  analyses 


Table  3.3 


Flight  Condition  Summary  for  Phase  II  Flight  Test 


MANEUVER  DESCRIPTION 

TAPE /FILE 

DURATION 

(SEC) 

ho 

(FT) 

M 

RANGE 

RANGE 

(DEG) 

8 

RANGE 

(DEG) 

Stall  with  Lateral  Doublets 

M01431/4 

58 

30K 

.3/ .5 

3/42 

-20/18 

Stall  with  Pitch  Doublets 

M01431/3 

74 

30K 

.35/. 45 

5/35 

-13/13 

Stall  with  Rudder  Doublets 

M01431/5 

87 

30K 

.32/. 36 

6/42 

-22/15 

Pitch  Doublets 

M01431/6 

44 

37K 

.4/. 49 

2/40 

-2/3 

Aileron  Doublets 

M00599/3 

42 

30K 

.4/. 45 

8/25 

-5/7 

Rudder  Doublets 

M00599/4 

57 

30K 

.49/. 52 

-8/25 

-20/25 

1G  Sequential  Doublets 

M00S99/5 

53 

3  OK 

.45/. 60 

0/33 

-15/10 

and  stall  and  departure  maneuvers  with  three-axis  doublets.  Two  flights  were 
flown  in  the  cleaning  loading  in  configurations  Cruise  (CR),  Take-off  (TO), 
and  Power  Approach  (PA)  to  update  thn  existing  clean-loading  simulation.  One 
flight  was  flown  in  a  high  stability  index  loading  and  the  same  configurations 
to  extend  the  simulation  capability  to  include  store  loadings.  This  unit 
stability  index  number  has  been  assigned  to  a  particular  weapon  or  piece  of 
suspension  equipment  as  a  measure  of  the  destabilizing  effect  of  that  item. 

The  high  stability  index  tests  were  performed  to  give  quantitative  estimates 
of  the  performance,  stability,  and  control  of  the  A/C  in  that  particular 
loading.  Table  3.4  presents  a  summary  of  the  Phase  III  test  conditions. 
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Table  3.4 

Phase  III  Flight  Test  Conditions 

> 

_  • 

FLIGHT  RECORD  CONFIGURATION 


TO 

TO 

TO 

PA 

PA 

PA 

PA 

CR  (M  -  .7) 

CR  (M  «  .8) 

CR  (M  -  .9) 

HSI  TO 
HSI  TO 
HSI  TO 
HSI  PA 
HSI  PA 
HSI  PA 

HSI  CR  (M  -  .7) 
HSI  CR  (M  -  .8) 


MANEUVER 


LONG.  SINE  LOW  AOA 

LONG.  SIN  MED  AOA 

RUDDER  SINE  MED  AOA 

LONG.  SINE  LOW  AOA 

RUDOER  SINE  LOW  AOA 

LONG.  SINE  HIGH  AOA 

RUDOER  SINE  HIGH  AOA 

LONG.  +  AILERON  +  RUDDER  SINE  LOW  AOA 

LONG.  SINE  LOW  AOA 

LONG.  SINE  LOW  AOA 

LONG.  SINE  LOW  AOA 

LONG.  SINE  MED  AOA 

LONG.  SINE  HIGH  AOA 

LONG.  SINE  LOW  AOA 

LONG.  SINE  MED  AOA 

LONG.  SINE  HIGH  AOA 

LONG.  SINE  LOW  AOA 

LONG.  SINE  LOW  AOA 


IV.  FLIGHT  TEST  DATA  ANALYSIS  AND  PREPROCESSING 


The  analysis  and  preprocessing  of  flight  test  data  is  a  major  element  of 
the  integrated  system  identification  procedure.  The  overall  objectives  of 
this  task  are  to  evaluate  the  suitability  of  the  flight  test  data  for 
identifying  aerodynamic  characteristics  and  to  develop  a  complete  set  of 
kinematically  consistent  measurements.  Completeness  is  required  to  ensure 
that  all  parameters  of  interest  can  be  identified.  Measurement  consistency 
has  a  significant  impact  on  parameter  identification  accuracy  since 
unaccounted-for  errors  will  bias  parameter  estimates. 

This  task  involves  five  steps: 

(1)  flight  data  evaluation; 

(2)  error  correction; 

(3)  data  filtering; 

(4)  data  reconstruction;  and 

(5)  instrument  calibration. 

These  steps  are  briefly  reviewed  below. 

4.1  FLIGHT  DATA  EVALUATION  AND  ERROR  CORRECTION 

The  flight  test  data  is  first  transferred  from  magnetic  tapes  to  binary 
disk  files.  For  each  flight  maneuver,  all  relevant  data  channels  are  plotted 
for  evaluation. 

Each  channel  is  inspected  for  sign  errors,  data  dropout,  and  wild 
points.  Sign  and  units  conversions  are  performed,  if  necessary.  Wild  points 
are  removed  using  an  interactive  computer  program  that  searches  for  adjacent 
data  points  differing  by  more  than  a  threshold  specified  by  the  user.  Data 
dropout  problems  can  be  solved  in  two  ways:  either  the  particular  variable  is 
reconstructed  from  other  measurements  or  if  the  data  loss  is  of  short 
duration,  the  maneuver  is  split  into  two  shorter  ones. 

The  second  part  of  the  flight  test  data  evaluation  concerns  the  quality 
and  completeness  of  the  measurement  set.  The  need  for  smoothing  the  data 
and/or  reconstructing  missing  or  noisy  measurements  must  be  assessed.  In 


addition,  the  frequency  content  of  the  maneuver  is  evaluated  to  ensure  that 
all  dynamic  modes  of  interest  are  properly  excited.  Based  on  this  evaluation 
several  different  maneuvers  (e.g.,  pitch,  roll,  and  yaw  sssp's)  are  combined 
to  improve  the  identification  accuracy. 

In  addition,  cross  plots  of  the  independent  variables  are  generated  to 
define  regions  where  the  nonlinear  models  can  be  identified.  Figure  4.1 
illustrates  an  independent  variable  cross  plot  which  has  been  generated  from 
Phase  II  flight  test  data.  This  crossplot  is  used  to  show  that  aerodynamic 
model  terms  which  are  a  function  of  a  and  e  should  be  identifiable  for  10 
<  a  <  30*  and  jgj  <10*. 

4.2  DATA  FILTERING 

A  Fast  Fourier  Transform  (FFT)  is  used  to  smooth  the  data  by  removing 
high-frequency  components.  The  FFT  algorithm  can  also  be  used  to  generate 
smooth  derivatives  for  any  desired  variable.  This  is  particularly  useful  to 

•  • 

generate  variables  such  as  P  or  i  which  are  not  usually  measured. 

To  be  effective,  the  filtering  process  should  remove  most  of  the 
high-frequency  noise  without  affecting  the  information  content  in  the 
frequency  range  of  interest.  Consequently,  the  choice  of  the  FFT  cutoff 
frequency  for  each  measurement  is  based  on  the  following  criteria: 

-  The  cutoff  frequency  must  be  at  least  four  times  greater  than  the 
highest  frequency  of  interest  (i.e.,  the  highest  dynamic  mode 
affected  by  that  particular  measurement). 

-  The  difference  between  the  original  and  filtered  signal  should 
appear  like  white  noise. 

Figure  4.2  illustrates  the  effect  of  different  cutoff  frequencies  on  the 
roll-rate  signal.  The  difference  between  the  original  roll-rate  signal  and 
the  filtered  signal  (with  a  cutoff  frequency  of  4  Hz)  is  shown  in  Figure  4.3, 
together  with  the  corresponding  curves  for  a,  8,  and  0.  The  proper  choice 
of  cutoff  frequency  is  even  more  critical  when  the  FFT  capability  to  generate 
derivatives  is  used  as  shown  in  Figure  4.4  (roll  acceleration  generated  from 
the  roll-rate  signal  of  Figure  4.2). 
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Figure  4.1 


Angle  of  Attack/Sideslip  Test  Data  Crossplot 
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gure  4.4  Effect  of 


4.3  DATA  RECONSTRUCTION 


Data  reconstruction  is  required  for  variables  that  are  not  measured  at 
all  or  not  measured  directly,  and  for  recorded  measurements  that  cannot  be 
used  due  to  instrumentation  or  data  acquisition  problems. 

For  the  F-4S  flight  test  data  analysis,  the  following  measurements  were 
reconstructed: 

(1)  The  roll  rate  signal  was  saturated  at  60  deg/sec  and  was 
reconstructed  from  the  roll  angle  derivatives  ($)  generated  by  the  FFT 
algorithm: 

P  ■  4  -  (Q  sin0  +  R  cos t)  tane 

(2)  The  angular  accelerations  were  generated  by  the  FFT  algorithm  from 
the  angular  rates.  The  equations  for  accelerometer  location  correction 
contain  terms  in  P,  Q,  and  R  and  these  derivatives  must  be  supplied  as  input 
when  the  corresponding  states  (P,  Q,  and  R)  are  not  propagated. 

(3)  Velocity  components  at  the  center  of  gravity  (U,  V,  and  W)  were 
reconstructed  from  the  nose-boom  measurements  (VT,  tt,  and  a).  For  example, 

VNB  -  VT  sine 


where 


£x,  lz  »  nose  boom  location  relative  to  the  c.g.  These  reconstructed 
velocities  are  used: 

(a)  To  obtain  the  initial  conditions  for  the  propagation  (i.e., 
integration)  of  u,  v,  and  w. 

(b)  To  enable  identification  of  only  a  subset  of  the  parameters. 
For  example,  to  identify  only  longitudinal  characteristics  (e.g., 

Cg(a)  and  CL(at ) ) ,  u  and  w  are  propagated  while  the  reconstructed 
values  are  used  for  v.  As  a  result,  CQ(a)  and  CL(a)  can  be 
identified  without  any  a  priori  knowledge  of  lateral/directional 
characteristics  such  as  C  . 

J  a 


(4)  Mass  and  inertia  characteristics  were  reconstructed  from 

measurements  of  the  fuel  quantity  in  each  of  the  aircraft's  eight  fuel  tanks. 

The  curves  relating  fuel  quantity  to  incremental  change  in  mass  and  inertia 

were  provided  by  the  NATC.  The  reconstructed  variables  were:  m,  I  I  , 

a  y 

V  *xz»  xcg»  Ycg*  and  zcg* 

4.4  INSTRUMENT  CALIBRATION 

The  instrumentation  system  errors  are  identified  so  that  sensor  signals 
will  be  kinematically  consistent.  Otherwise  bias,  scale  factors  and  time 
delays  in  the  sensor  signal  will  appear  as  errors  in  the  identified 
aerodynamic  parameters.  Kinematic  consistency  is  determined  by  comparing  the 
angular  positions  t,  e  and  4>  to  the  integrals  of  the  Euler  rates  which 
are  reconstructed  from  the  body  rates  P,  Q  and  R.  Similarly,  velocities 
(U,  V  and  W)  are  compared  to  the  integrals  of  the  accelerations.  Bias,  scale 
factors  and  time  delays  are  identified  so  as  to  reduce  the  comparison  errors. 

For  example,  pitch  rate-gyro  calibration: 


a©  =  e_ 


^1 

-/[  Q 


cos  -  R  sin  *m]dt 


%  -  bq >'<’  *  V 


(R.  -  br)/(l  *  Kr) 


where: 


■*m»  Ktn 


kq,  kr 
®mi  ^m 


true  pitch  and  yaw  rates 
measured  pitch  and  yaw  rates 
pitch  and  yaw  rate  bias  as 
pitch  and  yaw  rate  scale  factors 
measured  pitch  and  roll  attitudes 


j  and  k^  j  are  identified  by  the  SCIDNT  program  so  as  to 
minimize  ae.  Identification  of  the  time  delay  is  a  built-in  capability  of 
SCIDNT  and  is  done  by  shifting  the  Qm  input  until  a  e  achieves  a 


minimum. 


-  L 
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The  accelerometers,  angles  of  attack  and  sideslip  and  airspeed  sensors 
are  identified  in  a  similar  way  using  the  kinematic  equations  relating  Ax, 
Ay  and  Az  to  o,  b,  Vt  and  hp.  For  example,  vertical 
accelerometer  calibration: 

a  *  (a  -  b  )/(l  +  K  ) 
z..  v  zm  az"v  a z' 
aa  m 


zcg  zaa 


a,  MQ-PRH.  -(P+RQ)l„  *  (Q2  *  P2)i, 


«-/(«  .  Q  -  V  .  P  +  g  cose  cos«i  +  a  )  dt 


eg 


WNB  "  W  +  %  - 


Aa  ■  a„  -  tan 
m 


-i  /wnb\ 

\UNB  / 


baz  and  Kaz  are  identified  by  SCIDNT  so  as  to  minimize  Aa. 

It  should  be  noted  that  in  the  actual  process,  all  three  rate  gyros  and 
all  three  accelerometers  are  calibrated  simultaneously,  and  all  relevant 
measurements  are  used.  As  a  result  of  the  coupling  between  the  state 
variables,  additional  terms  (e.g.,  a  bias  on  the  vertical  gyro:  b^)  and 
cross  terms  (e.g.,  Kgp  and  Kazx)  can  be  identified. 

When  several  maneuvers  are  combined  together  and  processed 
simultaneously,  it  is  also  necessary  to  identify  maneuver-dependent  terms  for 
those  parameters  that  are  time  varying,  such  as  the  gyro  biases  (b^,  b@, 
and  b^)  and  rate-gyro  biases  (bp,  bq,  and  bR). 

The  complete  list  of  Instrument  calibration  terms  and  equations  is 
contained  in  Appendix  A,  and  Figure  4.5  presents  the  instrumentation 
calibration  model  which  was  Identified  from  Phase  II  flight  test  data. 


V.  SYSTEM  IDENTIFICATION  RESULTS 


This  section  presents  the  aerodynamic  and  propulsion  model  data  which 
were  identified  from  both  the  Phase  II  and  Phase  III  flight  test  programs. 
This  data  base,  along  with  the  mathematical  model  presented  in  Appendix  A, 
comprise  a  simulation  of  the  F-4S. 


5.1  LONGITUDINAL  CHARACTERISTICS  AND  EXTERNAL  STORE  EFFECTS 


5.1.1  General 


This  section  presents  the  F-4S  longitudinal  characteristics  at  low  to 
medium  angles  of  attack.  The  incremental  effect  of  external  stores  is  also 
shown.  Longitudinal  aerodynamic  models  were  identified  for  the  following 
configurations: 

-  Take-Off  (TO)  -  slats  and  landing  gear  extended,  trailing  edge 
flaps  deflected  30  deg,  and  inboard  leading  edge  flaps  deflected 
30”. 

-  Power  Approach  (PA)  -  slats  and  landing  gear  extended,  trailing 
edge  flaps  deflected  60*,  inboard  leading  edge  flaps  deflected 
30  . 

-  Cruise  (CR)  -  landing  gear  and  flaps  retracted,  slats  extend  and 
retract  automatically.  Two  CR  models  were  identified  for  Mach 
0.7  and  Mach  0.8. 

Two  models  were  estimated  for  each  configuration: 

-  A  clean  or  low  stability  index  (LSI)  loading. 

-  A  high  stability  index  (HSI)  loading  with  numerous  external 
stores. 

External  store  incremental  effect  was  derived  as  the  difference  between  the 
HSI  and  the  LSI  models. 


5.1.2  Instrumentation  Considerations 


Phase  III  flight  test  maneuvers  were  processed  for  this  phase.  Both  the 
LSI  and  HSI  flights  had  instrumentation  problems,  which  were  solved  as  follows. 

(1)  For  the  LSI  flights,  air  temperature  was  not  recorded,  and  a  cold 
standard  day  model  for  temperature  was  used.  This  impacts  the  air  density, 

Mach  number,  and  computations.  is  used  as  a  correction  factor 

for  the  thrust  model. 

(2)  For  the  HSI  flights,  the  engine  fuel  flow  and  flight  control  system 
(FCS)  gyro  data  were  not  recorded.  Since  fuel  flow  is  needed  in  the  gross 
thrust  model,  it  was  calculated  from  engine  RPM.  Figures  A. 3  and  A. 4  in 
Appendix  A  show  the  corrected  RPM  to  corrected  fuel  flow  curves  used  in  the 
analysis.  This  rleationship  is  only  valid  for  power  settings  less  than  mil 
power.  Since  the  FCS  vertical  and  rate  gyros  data  were  not  available, 
alternative  sources  of  these  measurements  were  explored.  The  measurements 
could  not  be  adequately  calibrated  by  using  instrumentation  pod  bank  angle, 
roll  rate,  or  pitch  rate.  However,  the  airdata,  noseboom,  and  inertial 
measurements  could  be  calibrated  when  roll  rate,  bank  angle,  yaw  rate, 
heading,  and  sideslip  were  assumed  zero  and  the  derivative  of  pitch  angle  used 
as  pitch  rate.  This  approach  is  valid  since  the  maneuvers  processed  for  the 
longitudinal  models  have  negligible  lateral-directional  excitation.  Table  5.1 
summarizes  the  sources  used  for  these  troublesome  measurements. 

The  instrumentation  problems  necessitated  identification  of  two  LSI  aero 
models.  The  first  and  more  accurate  model  is  based  on  the  FCS  and  fuel  flow 
data,  and  the  second  model,  LSI*,  is  based  on  the  instrumentation  pod  pitch 
attitude  and  calculated  fuel  flow.  The  LSI*  model  is  needed  to  determine 
store  effects  by  comparison  with  the  HSI  model.  The  LSI*  model  may  also  be 
used  in  the  future  if  fuel  flow  is  not  known. 


Table  5.1 


Measurement  Sources  for  LSI,  LSI*,  and  HSI  Analysis 


MEASUREMENT 

SOURCE  FOR  LSI  ANALYSIS 

SOURCE  FOR  LSI*  AND 

HSI  ANALYSIS 

ROLL  RATE 

FLIGHT  CONTROL  SYSTEM  (FCS) 
DATA 

ASSUMED  ZERO 

PITCH  RATE 

FCS 

SET  TO  4 

YAW  RATE 

INSTRUMENTATION  POD 

ASSUMED  ZERO 

BANK  ANGLE 

FCS 

ASSUMED  ZERO 

PITCH  ATTITUDE 

FCS 

INSTRUMENTATION  POD 

HEADING 

INSTRUMENTATION  POD 

ASSUMED  ZERO 

FUEL  FLOW 

RECORDED  FLIGHT  DATA 

CALCULATED  FROM  RPM 

AIR  TEMPERATURE 


ASSUMED  STANDARD  DAY  -15*C 


RECORDED  FLIGHT  DATA 


5.1.3  Results 


The  longitudinal  aerodynamic  characteristics  are  presented  as  follows. 

(1)  A  summary  of  the  data  processing  accuracy  is  presented  in  Figure 
5.1.  The  root  mean  squared  (rms)  errors  for  each  measurement  (Q,  a  ,  a,, 
aNB’  C|’  and  hp)  and  each  maneuver  are  shown. 

(2)  A  more  detailed  presentation  of  identified  model  accuracy  is  given 
in  Figures  5.2  to  5.4.  Each  figure  is  a  time  history  plot  for  one  of  the 
three  configurations  (TO,  PA,  and  CR).  Both  the  actual  recorded  test  data  and 
the  estimated  measurements  derived  from  the  identified  aerodynamic  models  are 
shown . 

(3)  The  identified  aerodynamic  coefficients  for  the  TO,  PA,  CR/M  =  0.7, 
and  CR/M  =  0.8  are  shown  in  Figures  5.5  through  5.10,  5.11  through  5.16,  5.17 
through  5.22,  and  5.23  through  5.28,  respectively.  Each  figure  presents  the 
variation  of  the  aerodynamic  coefficient  versus  wing  angle  of  attack  and  the 
corresponding  external  store  incremental  effect.  The  estimation  uncertainty 
is  shown  by  the  ±2a  bounds. 

The  lift  and  pitching  moment  coefficients  due  to  pitch  rate  are  shown  in 
Tables  5.2  and  5.3.  All  of  the  results  are  in  the  stability  axis  system  and 
are  referenced  to  c.g.  ■  33  percent  and  =  0. 
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Figure  5.1  Root  Mean  Square  Measured  Errors 
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Figure  5.2  Time  History  Plot  -  Takeoff  (HSI)  (Cont'd) 
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Figure  5.3  Time  History  Plot  -  Power  Approach  (HSI) 
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Figure  5.3  Time  History  Plot  -  Power  Approach  (HSI)  (Cont'd) 
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Figure  5.4  Time  History  Plot  -  Cruise  (HSI) 
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Figure  5.4  Time  History  Plot  -  Cruise  (HSI)  (Cont'd) 
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Figure  5.5  Lift  Coefficient  -  Takeoff 
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Figure  5.6  Drag  Coefficient  -  Takeoff 
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Figure  5.8  Pitching  Moment  Coefficient  -  Takeoff 
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Figure  5.12  Drag  Coefficient  -  Power  Aporoach 
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Figure  5.14  Pitching  Moment  Coefficient  -  Power  Approach 
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Figure  5.15  Pitching  Moment  Coefficient  due  to  Stabilator  - 
Power  Approach 
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Figure  5.18  Drag  Coefficient  -  Cruise  (M  -  .7) 


Table  5.2 


Pitching  Moment  Coefficient  Due  to  Pitch  Rate 


FLAP 

CONFIGURATION 

CLEAN  CONFIGURATION 

STORE  EFFECT 

C  *  2a 

mA 

q 

iC  *  2 a 
mq 

TO 

-6.40  ±  0.18 

-0.22  *  0.30 

PA 

-7.45  *  0.50 

-0.24  *  0.30 

CR  (M  -  0.7) 

-5.74  *  0.12 

-0.55  *  0.28 

CR  (M  =,  0.8) 

-6.28  *  0.08 

— 

Table  5.3 

Lift  Coefficient  Due  to  Pitch  Rate 


5.2  LONGITUDINAL  CHARACTERISTICS  AT  HIGH  ANGLE  OF  ATTACK 


This  sections  presents  the  F-4S  longitudinal  characteristics  at  low  speed 
(M  ■  0.4)  and  high  angle  of  attack  (5  to  35  deg).  The  test  maneuvers  used, 
taken  from  the  Phase  II  flight  test  are: 

-  Tape  3,  File  3  -  stall  with  pitch  doublet  at  30,000  ft  and  Mach 
0.4. 


-  Tape  3,  File  6  -  pitch  doublets  at  37,000  ft  and  Mach  0.4. 

The  results  are  presented  in  the  same  format  as  those  of  Section  5.3: 

(1)  Figure  5.29  shows  a  time  history  plot  of  one  test  maneuver 
with  both  the  actual  recorded  test  data  and  the  estimated 
measurement  derived  from  the  identified  aerodynamic  model. 

(2)  Figures  5.30  through  5.35  present  the  identified  aerodynamic 
coefficients  in  stability  axis,  referenced  to  c.g.  *  33 
percent  and  «s  =  0.  The  estimation  uncertainty  is  shown  by 
the  *2o  bounds. 
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Figure  5.29  (Continued) 
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Figure  5.29  (Continued) 
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5.3  LATERAL-0 IRECT IONAL  CHARACTERISTICS 


5.3.1  General 


The  F-4S  lateral /directional  characteristics,  presented  in  this  section, 
were  identified  in  two  distinct  phases.  The  initial  processing  covered  the 
full  low  to  high  angle-of -attack  range  and  was  based  on  the  Phase  II  flight 
tests.  The  final  processing  provided  more  accurate  estimates  for  the  low 
angle-of -attack  range  and  was  based  on  the  Phase  III  flight  tests  conducted. 
The  final  processing  was  based  on  improved  thrust  and  c.g.  location  models. 

5.3.2  Instrumentation  Considerations 


Initial  Processing 

In  the  initial  data  processing,  the  following  maneuvers  from  Flight  WD7 
were  used: 

-  Tape  2,  File  3  -  aileron  doublet 

-  Tape  2,  File  5  -  sequential  doublets 

-  Tape  3,  File  3  -  stall  and  pitch  doublet 

-  Tape  3,  File  4  -  stall  and  aileron  doublet 

-  Tape  3,  File  5  -  stall  and  rudder  doublet 

-  Tape  3,  File  6  -  pitch  doublet 

where  the  longitudinal  maneuvers  were  included  to  aid  i  instrumentation 
calibration.  All  the  maneuvers  were  flown  at  Mach  0.4  with  the  altitude 
varying  between  30,000  and  37,000  feet. 

A  net  thrust  propulsion  model  was  used  (see  Table  A. 6  in  Appendix  A  for 
details).  The  only  instrumentation  problem  was  with  the  roll-rate  gyro.  As  a 
result,  the  derivative  of  roll  attitude  was  used  to  reconstruct  roll  rate. 
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Final  Processing 

The  maneuvers  used  for  the  second  phase  were  the  pitch,  aileron,  and 
rudder  pumping  at  Mach  0.7  from  Flight  3,  Record  46.  A  more  accurate 
propulsion  model  representing  both  gross  thrust  and  inlet  ram  force  effect  was 
used  (see  Table  A. 7  in  Appendix  A).  An  updated  model  of  c.g.  location  as  a 
function  of  fuel  consumption  was  also  used.  The  maneuvers  covered  only  the 
low  angle-of -attack  range  and  a  linear  aerodynamic  model  (angle  of  attack 
perturbation  relative  to  o  =  5  deg)  was  sufficient.  Here,  too,  the 
derivative  of  roll  attitude  was  used  to  reconstruct  roll  rate.  A  time  history 
plot  of  the  second  phase  maneuvers  is  shown  in  Figure  5.36. 

5.3.3  Combined  Initial  and  Final  Estimates 

Before  the  initial  estimates  could  be  revised  using  the  final  processing 
results,  the  latter  had  to  be  corrected  for  Mach  number  effects.  The 
aerodynamic  coefficients  variations  with  Mach  number  between  Mach  0.7  and  Mach 
0.4  (at  o  ■  5  deg)  were  obtained  from  Ref.  43  (MDC  A2013).  The  initial 
estimate  at  a  =  5  deg  and  the  corrected  estimate  from  the  second  phase  (both 
at  Mach  0.4)  were  combined  according  to: 

C  -  K1C1  +  K2C2 


where 

C  =  revised  coefficient 

•  initial  coefficient  estimate  (Phase  II) 

C2  =  final  coefficient  estimate  (Phase  III)  corrected  for  Mach  effects 

o,  op  o2  =  the  corresponding  standard  deviation 


V  V.’ ■  ■■ 


The  revised  coefficient  (C)  was  then  used  instead  of  the  initial 
coefficient  as  the  o  ■  5  deg  knot  in  the  cubic  spline  representation  of  each 
aerodynamic  derivative. 

These  results  are  shown  in  Figures  5.37  through  5.44  in  the  following 
manner: 

(1)  The  top  plot  of  each  figure  shows  the  initial  estimate  as  a 
function  of  angle  of  attack  and  the  ±2o  uncertainty  bounds. 

(2)  The  estmate  from  the  second  phase,  at  ow  ■  5  deg,  are  also 

shown  on  the  top  plot.  The  square  (  □  )  is  the  value  for  Mach 

0.7  and  the  star  (*)  is  the  corrected  value  at  Mach  0.4. 

(3)  The  revised  estimate  is  shown  on  the  bottom  plot  of  each 

figure.  The  circle  (  O  )  at  aw  -  5  deg,  emphasizes  the  only 
spline  knot  that  was  changed.  The  other  knots  (at  15,  25,  and 
35  deg)  are  the  same  as  those  of  the  initial  estimate. 

NOTE:  The  lateral /directional  results  are  in  body  axis,  whereas  the 

longitudinal  results  are  in  stability  axis. 
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TIME  (sec) 

-  Estimate  from  Identified  Aerodynamic  Model 

oooo  Flight  Test  Data 


Figure  5.36  Time  History  Plots  -  Lateral/Directional 


TEST  INPUT:  AILERON  PUMPING 


(ft/sec  )  1S 


TIME  (sec) 
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Figure  5.36  (Continued) 
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TEST  INPUT:  AILERON  PUMPING 


Figure  5.36  (Continued) 


INITIAL  ESTIMATE 


Figure  5.37  Side  Force  Coefficient  due  to  Sideslin 


INITIAL  ESTIMATE 


Figure  5.45  Yawing  Moment  Coefficient  due 
to  Roll  Rate  (body  axis) 


—  Estimate 
.  .  Estimate  +2> 


5.4  PROPULSION  SYSTEM  MODELING  AND  IDENTIFICATION 


For  the  Phase  II  flight  test  processing,  a  simplified  model  of  the 
propulsion  system  was  used  to  separate  thrust  effects  from  the  applied  forces 
and  moments.  The  model  was  based  on  representing  the  engines  by  corrected  net 
thrust  as  a  function  of  corrected  fuel  flow.  This  model  was  identified  from 
flight  test  data.  The  result  of  this  identification  is  presented  in  Figure 
5.46.  The  difference  between  the  identified  and  preflight  models  for 
corrected  net  thrust  is  a  bias  in  fuel  flow  of  58  Ibs/hr  based  on  a  pressure 
altitude  of  hp  =  30,000  ft.  The  accuracy  of  the  flight  test  fuel  flow 
sensor  is  150  Ibs/hr. 

The  propulsion  model  was  revised  to  separately  represent  gross  thrust  and 
inlet  airflow  terms  for  the  data  processed  from  the  Phase  III  flight  test 
program.  These  models  were  not  identified.  Status  check  performance  data  for 
the  F-4S  engines  were  used  directly.  Details  of  this  model  are  presented  in 
Appendix  A. 

5.5  MODEL  VALIDATION 

The  validity  of  the  parameter  estimates  can  be  established  from  three 
different  considerations. 

(1)  Engineering  judgement:  are  the  estimates  reasonable  from  the 
point  of  view  of  general  agreement  with  preflight  predictions 

(2)  Estimation  uncertainty:  what  is  the  magnitude  of  the  +  2a 
bands  about  the  estimate 

(3)  Prediction  accuracy:  How  well  does  the  identified  model 
predict  flight  test  measurements  for  many  test  maneuvers  Does 
the  identified  model  predict  these  measurements  better  than  a 
model  based  on  preflight  parameters 

Figures  5.47  and  5.48  compare  identified  aerodynamic  characteristics  with 
preflight  predictions.  The  purpose  of  this  presentation  is  to  show  the 
reasonableness  of  the  identified  parameters.  The  identified  data  are  from  the 
Phase  II  analysis  and  repeat  data  presented  in  Section  5.3. 

Figure  5.47  presents  rudder  and  lateral  control  power  estimates  as 
functions  of  angle  of  attack.  The  flight  derived  control  power  estimates  show 


CORRECTED  NET  THRUST 


Figure  5.46  Identified  Net  Thrust  Model 


Figure  5.47  Comparison  Between  Flight-Identified  and  Preflight-Predicted 
Lateral/Directional  Control  Terms  (continued  on  next  page) 


good  agreement  with  preflight  data  in  terms  of  trends  with  angle  of  attack 
with  slight  differences  in  the  actual  value  of  the  derivative.  The 
uncertainty  of  the  lateral  control  derivatives  (C5  and  Cn  )  as  shown  by  the 

\ 

fanning  of  the  +2o  curves  is  due  to  the  lack  of  lateral  control  excitation  at 
low  and  high  o.  In  general,  when  the  uncertainty  boundaries  are  large  for 
some  parts  of  the  flight  regime,  this  information  can  be  used  for  planning 
additional  flight  test  conditions. 

Identified  values  for  roll  rate  and  yaw  rate  stability  derivatives  for 

the  F-4S  are  presented  in  Figure  5.48.  Good  agreement  is  evident  between 

preflight  prediction  and  identified  values  for  C  A  and  CoA  (for  a  < 

nP  P 

30).  The  identified  value  for  C  is  much  greater  than  the  preflight 

r 

prediction,  although  both  are  independent  of  angle  of  attack.  Significant 

differences  exist  between  identified  and  predicted  values  of  C0  for  a  > 

20*. The  mid-span  wing  pylon,  which  protrudes  beyond  the  wing  leading  edge  for 
the  flight  test  aircraft,  could  contribute  to  some  of  the  differences  between 
the  flight  identified  and  preflight  predicted  values  for  a  >  20*.  The  wing 
pylon  would  shed  a  vortex  that  would  alter  wing  aerodynamics  at  high  o. 

5.6  FLIGHT  MEASUREMENT  PREDICTIONS 

The  ability  of  the  identified  aerodynamic  and  net  thrust  models  to 
reproduce  flight  test  measurements  for  the  six  test  maneuvers  analyzed  from 
Phase  II  flight  test  data  is  tabulated  in  Figure  5.49.  Figure  5.50  summarizes 
the  root  mean  square  residuals  between  the  flight  test  measurements  and 
estimates  of  the  measurements.  The  predicted  measurements  are  generated  from 
a  simulation  of  the  F-4S  based  on  the  identified  aerodynamic,  thrust,  and 
instrumentation  models.  The  simulation  is  run  with  actual  flight  test  time 
histories  of  the  control  commands. 


RMS  RESIDUAL  FOR  COMBINED  SIX  MANEUVERS 
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00 
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• 
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q 
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2.22 

6.56 

% 

ft 
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n 

X 
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n 

y 
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nz 
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Note:  Measure  Residual  s  Actual  Measurement  -  Estimate  of  Measurement 


Figure  5.49  Effect  of  Model  Parameters  on  Measurement  Residuals 


-  PREDICTED 

oooo  MEASURED 


-  PREDICTED 

oooo  MEASURED 


To  gain  some  feeling  for  the  relative  worth  of  the  identified  model,  the 
F-4S  simulation  was  configured  with  the  preflight  aerodynamic  data  and  run 
with  the  flight  test  control  inputs.  The  resulting  measurement  residuals  are 
also  presented  in  Figure  5.49.  In  all  cases,  the  identified  model  provides  a 
better  explanation  of  the  flight  test  measurements.  Figure  5.50  compares  the 
time  histories  of  the  flight  test  measurements,  predictions  based  on  the 
identified  model  and  predictions  based  on  the  preflight  model  for  one  of  the 
flight  conditions  analyzed. 
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VI.  SUMMARY  AND  CONCLUSIONS 


The  results  presented  in  this  report  demonstrate  the  operational  status 
of  a  nonlinear  system  identification  <^ta  processing  technique  which  is 
currently  being  used  at  the  Naval  Air  Test  Center.  The  following  specific 
conclusions  are  noted. 


(1)  Nonlinear  system  identification  data  processing  techniques  can 
be  used  to  identify  aerodynamic,  propulsion  system  and 
instrumentation  calibration  models  from  a  common  set  of  flight 
test  conditions. 

(2)  Aircraft  response  predictions  are  improved  with  the  identified 
model  relative  to  prediction  based  on  a  preflight  model. 

(3)  Extraction  of  performance,  stability  and  control,  and  high 
angle  of  attack  characteristics  from  a  single  model  has  been 
illustrated. 

(4)  The  capability  for  identifying  nonlinear  aerodynamic  models  in 
a  format  compatible  with  preflight  predictions  has  been 
demonstrated. 

(5)  A  methodology  for  determining  tne  accuracy  of  the  parameter 
estimates  was  presented. 

(6)  Dynamic  test  techniques,  which  require  nonlinear  system 
identification  data  processing  techniques,  can  improve  test 
productivity. 

(7)  It  is  feasible  to  identify  both  aerodynamic  and  propulsion 
model  data  from  a  common  set  of  flight  test  data.  Only  a  very 
simplified  propulsion  model  was  Identified  in  the  analysis 
covered  by  this  report.  To  aid  the  identification  of 
propulsion  models,  it  is  recommended  that  the  test  conditions 
include  throttle  transients  and  that  measurements  of  fan  inlet 
pressure  (Pj2),  turbine  exit  on  tail  pipe  pressure  (PT56)* 

and  nozzle  area  (Ag)  also  be  included.  As  could  possibly 
be  reconstructed  from  nozzle  actuator  displacements. 
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APPENDIX  A 

F-4S  MATHEMATICAL  MODEL 


This  appendix  describes  the  F-4S  mathematical  simulation  model 
implemented  by  the  SCIDNT  system  identification  computer  program.  The 
simulation  models  the  following  aspects  of  aircraft  dynamics: 

-  rigid  body  nonlinear  three  degree  of  freedom  rotation  and  three 
degree  of  freedom  translation; 

-  aerodynamic  and  propulsion  system  forces  and  moments; 

-  winds  (north,  east,  and  down  components);  and 

The  simulation  also  models  sensors  which  are  commonly  employed  during 
aircraft  flight  testing.  These  sensors  are 

-  angle  of  attack  and  sideslip  vanes  mounted  on  a  nose  boom; 

-  impact  pressure; 

-  true  airspeed; 

-  pressure  altitude; 

-  three  orthogonal  axis  accelerometers; 

-  three  orthogonal  axis  angular  rate  gyros; 

-  vertical  and  directional  gyros. 

Provision  is  also  made  for  modeling  scale  factor,  misalignment,  and  bias 
errors  in  all  sensors. 

The  following  is  a  list  of  the  subsections  describing  the  simulation 
model : 


A. 1  STATE  EQUATIONS 
A. 1.1  Definitions 

A. 1.2  Fuselage  Rigid  Body  Kinematics  and  Dynamics 
A. 1.3  Applied  Forces  and  Moments 
A. 1.3.1  Modeling  Approach 
A. 1.3. 2  Expansion  Variables 
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A.2  Measurement  Equations 
A. 2.1  Inertial  Sensor 
A.2. 2  Air  Oata  Sensors 
A. 2. 3  Sensor  Moment  Arms 

A.l  STATE  EQUATIONS 

A.  1.1  Definition  of  States 

Table  A.l  lists  the  10  states  used  in  the  simulation  model.  This  list  is 
made  up  of  three  basic  classes  of  states 

-  fuselage  translational  and  rotational  motion, 

-  attitude  dynamics,  and 

-  altitude  dynamics. 

The  units  listed  are  in  generic  length  (L),  mass  (M),  time  (T),  and  angle 
(rad).  The  program  works  effectively  for  any  consistent  set  of  units  such  as 
International  or  English.  The  geometric  quantities  cgx,  cg2,  XREF, 
and  ZREF  will  be  scaled  by  1/12  before  being  used  to  compute  Xcg  and 
Z__  (see  A.l. 3.1).  This  will  simplify  the  use  of  inches  for  the  first  four 
terms  and  feet  for  the  second  two.  Dimensional  constants  such  as  air  density 
must  also  use  consistent  units.  Table  A.2  illustrates  these  definitions. 


A. 1.2  Fuselage  Rigid  Kinematics  and  Dynamics 


The  detailed  differential  equations  which  describe  the  rigid  body  states 
and  the  kinematic  relationships  are 

-  Body  Fixed  Inertial  Velocities 
$  ■  W-“Q  *  g  TU.3)  *  <±>  Fx 

§  -  WP-UR  ♦  9  Ttt,3>  *  (£)  Fy 
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Table  A.l 
State  Definitions 


DESCRIPTION  UNITS 


x  body  axis  component  of  inertial  space  ft/sec 

velocity  of  body-fixed  origin 

y  body  axis  component  of  inertial  space  ft/sec 

velocity  of  body-fixed  origin 

z  body  axis  component  of  inertial  space  ft/ sec 

velocity  of  body-fixed  origin 

x  body  axis  component  of  angular  rotation  rad/sec 

(roll  rate) 

y  body  axis  component  of  angular  rotation  rad/sec 

(pitch  rate) 

x  body  axis  component  of  angular  rotation  rad/ sec 

(yaw  rate) 

Euler  roll  angle  rad 

Euler  pitch  angle  rad 

Euler  yaw  angle  rad 

Altitude  (positive  up)  ft 


Table  A. 2 

Dimensional  Parameters 
CONSTANTS 


SYMBOL 

DESCRIPTION 

UNITS 

9 

Acceleration  due  to  gravity 

fps2 

p 

Density  of  air 

slug-ft3 

m 

Mass  of  aircraft 

slugs 

lx 

Roll  moment  of  inertia  =  (y2+z2)dm 

slug-ft2 

h 

Pitch  moment  of  inertia  =  (x2+z2)dm 

slug-ft2 

h 

Yaw  moment  of  inertia  =  (x2+y^)dm 

slug-ft2 

*xy 

x-y  product  of  inertia  =  xy  dm 

slug-ft2 

Ixz 

x-z  product  of  inertia  =  xz  dm 

slug-ft2 

y-z  product  of  inertia  «  yz  dm 

slug-ft2 
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S  -  UQ-VP  *  9  T  (3,3)  *  (i)  FZ 

LCD 


Body  Axis  Rotation  Rates  (I  =1*0) 

xy  yz 

n  ■  C1RQ  *  C2P«  *  C3"rx  *  C4hTZ 


§  -  CgRP  *  C6(R2-P2)  *  C7MTy 
37  *  C8Pt|  *  C9RQ  *  C4MTX  *  C10MTZ 


-  Euler  Angles 

^  »  P  +  Q  sin$  tan©  +  R  cosd  tan© 
«  Q  cos$  -  R  sintf 
^  ■  (Q  sin$  +  R  cos4)/cos© 


-  Altitude 
dh 

Tfr  *  -T  (1,3)  U  -T  (2,3)  V  -T  (3,3)  W 
ax  L2B  L2B  L2B 

In  these  equations,  the  matrix  T^g  is  a  direction  cosine  matrix 
which  relates  the  local  earth-fixed  coordinate  system  to  the  fuselage 
body-fixed  coordinate  system.  The  gravity  vector  is  aligned  with  the  z  axis 
in  the  earth-fixed  frame.  The  x  y  z  axes  of  the  body-fixed  frame  are 
aligned  with  the  reference  geometric  axes  of  the  fuselage.  The  nine  elements 
of  TL2b  are: 

T  (1,1)  .  cos©  cosi> 

L2B 


T  (2,1)  »  sin0  sin©  cos^  -  cosd  sin^ 


T  (3,1)  -  sind  si  nip  +  cosd  sin©  cosip 
L2B 

T  (1,2)  *  cos©  simp 
L2B 

T  (2,2)  «  cosd  cosip  +  sind  sin©  simp 
L2B 

T  (3,2)  =  cosd  sin©  simp  -  sind  cosip 
L2B 

T  (1,3)  -  -sind 
L2B 

T  (2,3)  =  cos©  sind 
L2B 

T  (3,3)  =  cos©  cosd 
L2B 

For  any  vector  x,  the  components  in  the  local  frame  and  the  components  in 
the  body  frame  are  related  by 

•^BODY  "  £TL2B^  £XWal 
AXIS  AXIS 


The  inertia  terms  C,,  i»l  to  10  are 
r  f'v-'z)  ‘zV 

A  T  T  T  t 


‘x‘z  -  ‘xz 


r  (Ix-!,*'z>  ‘xz 

c*  ■  ,  nr  i 

X  z  xz 


I. 


3  I  I  -I 2 
x  z  xz 


I 


xz 


4  I  I  -I  2 
X  z  xz 


M  » 
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I 


XZ 


y 

1.0 


.  <V'y>  Vxz' 

8  ii7 — 

x  z  xz 

("y-'z-'xl  'xz 
C  .  —  —  —  — 2 — 

9  ‘x'z-'xz 


1. 


'10  “  7 

iU  I  I  -I  c 
X  z  xz 


The  use  of  inertia  terms  in  this  form  assumes  symmetry  through  a  vertical 

plane  of  the  aircraft  (i.e.  IYU*IU,»0) .  The  aircraft  inertial  tensor 

*y  yz 

terms  are 

Ix  roll  moment  of  inertia  =  /(y2  +  z2)dm 
Iy  pitch  moment  of  inertia  =/  (x?  +  z2)dm 

I2  yaw  moment  of  inertia  =  /(x^  +  y2)dm 

IXz  x-z  product  of  inertia  =  /xz  dm, 
and  are  defined  for  the  body  axis  system. 

A. 1.3  Applied  Forces  and  Moments 

A. 1.3.1  Modeling  Approach 

The  applied  force  and  moment  terms  are  formulated  from  aerodynamic  and 
propulsion  system  contributions  (see  Table  A. 3).  To  permit  multimaneuver 
identification  processing,  the  moment  models  are  referenced  to  a 


t 
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Table  A. 3 


Applied  Forces  and  Moments  (Body  Axis) 


Axial  Force 


FAx  +  FTx 


Side  Force 

Fy  "  FAy  +  FTy 
Vertical  Force 


Fz  -  FAz  *  FTz 


Rolling  Moment 


M  *  M.  +  Mt  -  F  .  Z 
x  Ax  Tx  y  eg 


Pitching  Moment 

My  ■  %  *  MTy  -  Fy  •  Zcg 
Yawing  Moment 


F_  •  Xrn 

z  eg 


M  =  M  +  Mt  -  F  .X 
z  Az  Tz  y  eg 


FC)*  MC) 
FA(*)’  %•) 

T( *) '  MT( *) 

X  ,  Y  ,  Z 
ca  ca’  ca 


=  total  applied  forces  and  moments 
=  aerodynamic  forces  and  moments 
=  propulsion  forces  and  moments 
=  moment  arm  transfers  between  the  MRC  and  the  eg 


user-specified  "moment  reference  center"  and  then  are  corrected  by 
conventional  moment  arm  transfers  for  the  c.g.  location  of  each  maneuver. 

The  moment  arm  transfers  are  defined  by 

*cg  “  (cgx-Xj^p)  (^) 

*cg  *  -0833  (cg^-Z^pp)  (ft) 

where  cgx:  eg  Body  Station  location  (inches) 

cg^:  Water  Line  location  (inches) 

XREF:  Moment  Reference  Center  Body  Station  (inches)  - 
longitudinal 

YREF  ■  0 

Zref:  Moment  Reference  Center  Water  Level  (inches)  - 
vertical. 

The  aerodynamic  force  and  moment  models  are  defined  in  terms  of  nondimensional 
coefficients  which  are  dimensional ized  by  the  appropriate  reference  geometry 
and  by  the  dynamic  pressure.  The  aerodynamic  nondimensional  coefficients  are 
further  expanded  by  nondimensional  derivatives  represented  by  single  cubic 
spline  In  a  and  a.  The  aerodynamic  derivatives  are  defined  in  Table  A. 4 
for  the  body  axis  system  and  in  Table  A. 5  for  the  stability  axis  system. 

The  propulsion  forces  and  moments  are  defined  in  terms  of  corrected  fuel 
flow  and  corrected  engine  RPM.  Two  models  were  used  in  the  F-4S  analysis. 

The  first  was  a  simplified  model  using  net  thrust  as  a  function  of  corrected 
fuel  flow.  The  second  was  a  more  accurate  model  using  gross  thrust  and  inlet 
ram  force  as  a  function  of  corrected  fuel  flow  and  corrected  engine  RPM. 

These  two  models  are  defined  in  Tables  A. 6  and  A. 7,  respectively. 
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Table  A. 4 


Aerodynamic  Model  in  Body  Axis 


qS  Cv 


qs  C, 


qS  C, 


qS  b  C. 


REF 


qS  c  C 


m, 


REF 


qS  b  C 


'REF 


Cx<“>  *Cx  ,<“>  •  I ‘si 


C  («)  .  B  +  •  C  (a)  •  +  C  (a)  .  5 


X,  ’ 

\  ^ 


R 


+  C  .  R 
yR 


C,(a)  ♦  Kr  .  C  (a)  .  «c  +  C,  .  q 


Table  A. 4  (Continued) 


Sef  ‘  Va)  • '  *  V  •  Sr  *  *  C\M  ■  p*  *  \M 


c»<«>  *  Cm >>  •  ‘s  *  cm<?>  •  I'l  *  Cm  •  «*  *  V“!  '  '•*' 


l<R 


Cn<-  9>  *  Cn.<“>  •  SR  *  Cn  (•>  •  4A  *  C„i»>  •  **  +  V  •  R* 


Notation: 


a  single  cubic  spline  in  a 


Z  C(o  *  a. )  f . (a) 

i-1  1  1 


C(a,  b)  »  a  bicubic  spline  in  a  and  b,  which  is  antisymmetric  in  b 


k  k 


Z  Z  C(a  =»  Oj  ,  B  *=  bJ  fi(a)  g.(s) 
i-1  j-2  1  J  i  J 


-  C(a,  -B) ;  C(a,  0)  -  0  (for  g  <  0) 


(for  b  >  0) 


Table  A. 5 

Aerodynamic  Model  in  Stability  Axis 


FSx  =  ^SwCD 
FSy  =  ^SwCy 
FSz  =  "  ^SwCL 


k 


f 


M_  =  qS  b  C. 
Sx  Hwwi 


REF 


Mc  =  qS  c  C 
Sy  w  w  mREF 


Mc  =  qS  b  C 
Sz  H  w  w  n 


REF 


where 


FS(.),  Ms(.)  =  aerodynamic  forces  and  moments  in  stability  axis. 


As  the  equations  of  motion  are  solved  in  body  axis,  the  following 
stability  axis  to  body  axis  transformation  is  used: 


FAx  *  FSx  cos  aB  ”  FSz  sin  °B 


*  '« 


FAy  *  FSy 


FAz  31  FSx  sin  °B  +  FSz  cos  “B 


Ax  =  "Sx 


B  ‘  Sz 


MAy  “  MSy 


M.,  =  Me v  Sin  Op  +  Me,  COS  aR 


Table  A. 5  (Continued) 


The  aerodynamic  coefficients  are  expanded  as  follows: 


CD  “  CD^a^  +  CDi  W  *  ^S1 

usi 

S  ■  CyW  •  *  *  '  •„  »  Cy|a)  •  l.Al  *  ■  P  *  C,r  •  S 

CL  .  CL(.)  *  KC^  •  Cm(.)  *  *s  *  •  Q 

Cl  =  CA  (a)  *  e  +  C£(a)  •  5r  +  C£(o)  *  SA  +  C£(a)  *  P 
REF  8  fiR  6^  P 

+  C4(a)  *  R 
r 

CmR£F  =  Cm(a)  +  Cm^  *  6S  +  Cm(“>  *  * 

C"REF  ’  '"a1"1  ' 8  *  C"iRa>  •  s«  *  M;1 ' *  c4a) ' ' 

+  c  (a)  •  R 
"r 


NOTE:  Though  physically  and  numerically  different,  the  same  symbols  are  used 
for  body  axis  and  stability  axis  lateral  rolling  and  yawing  moment 
coefficients  (C^p  and  CngEp).  The  side  force  coefficient 
(Cy)  is  physically  and  numerically  the  same  in  both  axis  systems. 


Propulsion  Model  Using  Net  Thrust 


where 


and 


FTx  -  <fNL  *  \>  cos  ®T 
FTy-0 

FTz  -  -  (FNl  *  V  si"  ®T 
MTx  '  |lNy!  \  *  FNr>  s1"  sT 
%  '  *NZ  (F»l  *  \)  cos  ®T 
MTz  ■  l*Hyl  (FNl  *  F«R>  cos  »T 

Fn  ,  fN  *  left  and  right  engine  net  thrust 
\  "a 

«.N  ,  SLf,  ■  lateral  and  vertical  net  thrust  moment  arm 
y  nz  relative  to  MRC 

©t  «  thrust  line  inclination  relative  to  x-axis 


where 

fn  ■  corrected  net  thrust  versus  corrected  fuel  flow  curve 
identified  from  flight  test  data  (see  Section  5.4). 


Table  A. 7 

Propulsion  Model  Using  Gross  Thrust 


FTx  “  i«L,R  (FG.  cos  ®T  "  FI.  cos  “B  cos  B> 
FTy  "  "  i«L,R  (FI.  Sin  b) 


FTz  ■  -  2  (Fg  sin  ©T  +  Fj  sin  ag  cos  a) 


£ 

MTx  "  i-L,R  £FG.  *G  sin  ®T  +  FI.(£I  sin  B  ~  ll  cos  B  sin  aB^ 


Mjy  *  j £Fq ^ (  2.q  sin  0j  ^  2-q  cos  8j) ^  F j  ( sin  <xg  cos  8 
1  yi  zi  1  xi 


-  Hr  COS  a„  cos  b)] 
2i 


MT_  =  i  f  o  l-  Fr  lc  cos  ®T  +  Fi  (-  sin  e  +  tT  cos  aR  COS  b)] 

TZ  ,.L,R  G,  Gyj  T  I,  IXi  Iy 


where 


V  \ 

F.  ,  F, 
lL  lr 

y(  )  z(  i 


*i  £i 

x(  )  y(  )  z(  ) 


=  left  and  right  engine  gross  thrust 


left  and  right  inlet  ram  force 

gross  thrust  moment  arms  relative  to  MRC 


inlet  ram  force  moment  arms  relative  to  MRC 
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Table  A.7  (Concluded) 


corrected  air  mass  flow  versus  corrected  RPM  curve 
supplied  by  the  manufacturer  (see  Figure  A. 2) 


For  some  of  the  flights,  the  fuel  flow  measurement  was  not  available. 
Calculated  fuel  flow  as  a  function  of  corrected  RPM  was  used  instead 
(see  Figures  A. 3  and  A.4). 


A. 1.4. 2  Expansion  Variables 


The  expansion  variables  are  the  Independent  arguments  for  the  applied 
force  and  moment  terms.  These  variables  fail  into  four  classes  and  they  can 
be  calculated  as  a  part  of  the  simulation  or  read-in  as  "control  variables." 

-  Airspeed,  Angle  of  Attack,  Sideslip 

-  Body  Rotational  Rates 

-  Control  Positions 

-  Atmospheric  Variables 

-  Airspeed,  Angle  of  Attack  and  Sideslip 

The  airspeed  components  are  reconstructed  from  the  inertial  velocity 
components  (U,V,W)  and  the  body  axis  wind  velocity  components  (l^,  vw, 

WW>: 

UA  “  U  -  UW 


wA  -  w  -  ww 

The  body  axis  wind  velocity  components  are  derived  from  wind  components 
which  are  defined  in  an  earth  fixed  axis  system  (i.e.,  VN  (north),  VE 
(east),  VQ  (down)).  The  transformation  from  earth  to  body  axis  wind 
components  is 


(u”l 

1 

/  y 

‘  -  CTL2B^  ' 

* 

(ww  J 

1 

v 

Airspeed:  Vy  -  /UjJ  +  VA  +  WA 


angle  of  attack: 


angle  of  sideslip: 
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